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Summary 


A model is developed for the formation and propagation through 
the lower corona of the loop-like coronal transients in which mass 
is ejected from near the solar surface to the outer corona. I'he mass 
ejections are simulated with numerical solutions of the time-dependent, 
two-dimensional, dissipationless, magnetohydrodynamic equatlon.s of 
motion. The objective is to obtain a better understanding of the 
physical processes responsible for and occurring in the mass ejections. 

The main difference between this model and previous similar 
models is that we assume that the initial state for the transient is 
a coronal streamer. In previous models the initial state has been 
taken as a hydrostatic state with a force-free magnetic field. In the 
coronal streamer, the atmosphere is not stationary and the magnetic 
field is not frrce-free. A coronal streamer consists of closed magnetic field 
lines near the solar surface with overlying’ and adjacent open field lines. The 
plasma is stationary in the closed region and flows outward in the open region. 

Coronal streamers are observed to last for long periods on the 
solar surface and, the\|efore, must represent a quasi steady~s.*ate 
solution to the complete time-dependent equations. Since a self- 
consistent analytic solution for coronal streamers does not exist, 
a numerical solution is constructed. The initial state for the 
streamer is a polytropic, hydrodynamic solution to the steady-state 
radial equation of motion coupled with a force-free dipole magnetic 
field. The numerical solution of the complete tljne-dependent equations 
then gradually approaches a stationary coronal-streamer configuration, 

The streamer configuration becomes the initial state for the coronal 


transient, The streamer and transient simuXatlonB are performed 
completely Independent of each other. The transient is created by 
a sudden increase in the pressure at the base of Che closed-field 
region in the streamer configuration. Both coronal streamers and 
coronal transients are calculated for values of the plasma beta 
(the ratio of thermal to magnetic pressure) varying from 0.1 to 100, 

We compare our results with similar results obtained using the 
previously-developed model mentioned above in which the initial 
atmosphere is stationary and the magnetic field is force-free. We 
also compare both our results and the previous results with the 
observed characteristics of loop-lJke, mass ejection coronal transients 
in order to evaluate the relative merits of the two models. 
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I. INTRODUCTION 


Coronal mans ejections are transients or disturbances in the corona 
in which material appears to escape from the Sun. The data obtained 
during the Skylab mission, by both on-board and ground-based instrumenta- 
tion, have contributed significantly to the understanding of coronal 
mass ejections. For several ejections the coverage was complete enough 
that observations of the coronal transient (seen in white light and radio) 
could be associated with observations related to the solar phenomena, 
usually eruptive prominences or flares (seen in XUV, EUV, and X-ray), 
responsible for the transient (Hlldner et al. , 1975a, bj Rust and Hildner, 

1976; Dulk et al. , 1976; Poland and Munro, 1976; Schmahl and Hildner, 

1977). This complete data analysis and establishment of a cause-effect 

relationship are crucial in determining the physical processes in mass 

ejections. These associations are particularly informative for the 

purpose of developing dynamic models when the data can be used to infer 

the thermodynamic properties (temperature, density and pressure), 

velocities and magnetic fields of at least some portions or features 

of either the solar event or the transient. If the simulation can be 

made to reasonably reproduce the observed physical properties, then 

since all the above properties are included in the simulation, information will 

be obtained on the physics of the solar event and transient for which 

observations are not available. 

I 

Several models have been developed to simulate the dynamics of 
the transient (see, e.g., Nakagawa et al. , 1975; Steinolfson and 
Nakagawa, 1976; Steinolfson et al. , 1978; Wu et al. , 1978). These 
models do clarify some basic physics of the observed features of mass 


ejecciona*, however, they ore not able (by the very nature o£ the modeXa) 
to simulate some other Important features. The shortcomings of the 
current models (which will be discussed in Section V) are such 
that it is questionable whether they are capable of simulating 
the detailed structure of mass ejectlonsj even though they may be 
applicable to other types of atmospheric transients. In this report 
we present a model which lead.s to a more realistic simulation of the 
observed features of mass ejections and, consequently, to a better 
understanding of the related physical processes. I’he particular 
class of ejections which we simulate is the loop-llke mass ejections 
that tend to be associated with eruptive prominences (Hildner, 1977). 

The equations used in the model and the procedure followed in 
solving the equations are discussed in Section II. The complete 
simluatlon consists of two parts. In the first part a coronal- or 
helmet-streamer configuration is established; this is discussed in 
Section III. A coronal-streamer configuration consists of closed 
magnetic loops on the solar surface that lie beneath and adjacent to 
an open-field region. In the second part a perturbation in the thermal 
pressure is Introduced at the base of the closed loops which simulates 
the solar event. The numerical solution of the time-dependent equations 
then simulates the resulting coronal transient. The coronal transient 
is discussed in Section IV. The model is compared with a previous model 
in Section V and with observations in Section VI. The results are 
summarized and the conclusions that can be made from this study are 
discussed in Section VII, 


11. EQUATIONS AND SOLUTION PROCEDURE 

The solar atmosphere Is assumed to behave as a single fluid with 
negligible dissipative effects. With these assumptions the time-dependent 
magnetohydrodynamic (MHD) equations that describe atmospheric flows in 
the meridional plane can be written in MKS units as follows! 
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where the dependent variables are the density p, radial velocity u, 
meridional velocity v» pressure p, radial magnetic field and 
meridional magnetic field Bg. The independent variables are the radius 
r and the latitude 0. The constants are the poly tropic index Y» 
magnetic permeability y, solar gravitational constant G, and solar 
mass Mg. 
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fhe region in which the solution to the above equations 1§ desired 
is shown in Figure 1; l.e, , the region bounded by the solar surface and 
5 solar radii (R^) In radial distance and by the equator and the pole In 
meridional distance, the solution Is assumed fcy be. symmetric about the 
equator. The equations are solved numerically using a modified Lax- 
Wendroff difference scheme given by Rubin and Burst^in (1967) . The grid 
spaclngs used are Ar ** O.l R^ and AD“2.5 deg, The time step la chosen 
to be the maximum allowable from the usual stability criterion for 
Eulerlan difference schemes} l.e,, At »• min (At^, At^) where Atj, *» Ar/jXrj 
and Atg rAO/IXg] and |X^| (jXgi) is the maximum eigenvalue (the sum of 
the fluid velocity and the characteristic velocity) In the radial (meridional) 
direction. A smoothing term suggested by Lapldus (1967) Is used to re- 
duce numerical oscillations. The amount of the smoothing is controlled 
by a constant, Lapiduii vsed a value of 4; we use a value of 2. The 
initial conditions and the boundary conditions for this initial-boundary 
value problem for the coronal streamer and the coronal transient calcu- 
lations are dlepussed in Sections III and IV, respectively, 

Wliert Equations (1) arc written In conservation form, as they must be 
in order to solve ivumerlcally , they contain singularities at the pole due to 
terms of the form v/sln 0 and Bg/sln 0. Both the numerator and denominator of 
these relations vanish at the pole, In the numerical solution these singu- 
larities were overcome by setting the values of the above two relations at the 
pole (point 2 In Figure 1) equal to their respective values at the grid point 
next to the pole (point 3) . 

The above equations can be non-dimetislonalized with respect to some 
selected reference values in such a way that the only parameters are 
Y,,.13 (the ratio of thermal to magnetic pressure), and a parameter 
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involving the reference temperature. We aeleoted the reference valuet* 
ae the Initial values at the solar surface at the equator. Since > nml 
the reference temperature are kept constant, the only parameter remaining 
in the equations is A number of theoretical studies have demonstrated 
the Importance of the value of 0 (see, e.g., Steinolfson and Dryer, 1978; 
Nakagawn et al. , 1976; Steinolfson ot al,, 1978). The importance of 0 
in the present problem is investigated by examining eoronnl -streamer nnd 
coronal-transient structures for 0 varying from 0.1 to 100, The results 
for 0 « 4 and 0 ■ 0,5 are discussed in the most detail. 

A listing of the computer code used to numerically solve Equatlona (1) 
is included In Appendix A along with a description of the input variables 
for the code, A separate code is used to produce the computer-gansmtod 
plots presented in this report, and this computer code and o description 
of the input variables for it arc given in Appendix B. Appendix C 
contains a listing of a sample run for a coronal-streamer simulation. 

All of the computer simulations were performed on the CRAY-X computer at 
the National Center for Atmospheric Research in Boulder, Colorado, 
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III. CORONAL STREAMERS 

A schematic of the magnetic field lines in a coronal streamer is 
shown in Figure 2. The fluid velocity is essentially zero inside the 
cloBod-flald region with 11 the outward flow occurring in the open-field 
region. A current sheet exists on the equator in the open region where 
the magnetic field reverses direction. Coronal streamers are long-lived 
structures (often lasting for days) and, hence, r'MSt represent, in some 
sense at least, a configuration approaching a solution to the steady 
state form of Equations (1). However, such « steady-state solution would 
be extremely difficult to obtain. Pneuman and Kopp (1971) have con- 
structed a steady state solution, but they assumed that the field lines 
were known a priori and neglected the interaction between the field 
lines and the fluid. Our approach is to start with an initial state 
as close to a coronal streamer as feasible, numerically solve the time- 
dependent equations and let the solution relax with time to a coronal- 
streamer configuration. This approach v/as used by Endler (1971) and 
Weber (1978) for a constant temperature medium} an assumption which 
we do not make. The final coronal-streamer configuration serves as 
the initial state for the coronal transient discussed in Section IV. 

A. Initial State and Boundary Conditions 

The thermodynamic variables and the velocity arc given initially 
by a radial, hydrodynamic solution to the time-independent form of Equations 
(1); i.e., a Parker-type solution (Parker, 1963). The reference values 
for the thermodynamic variables in the initial state at 1 R^^ are 
temperature T »1.8 x 10^ K and electron number density n “ 2.25 x 10^ cin~^, 
where p *> nnip and is the proton mass. With these values the only 
solutions to the steady-state hydrodynamic equations with the critical 


point beyond 1 are obtained for 1 < Y < 1.1. The larger the value of 
Y. the more rapidly 3 increases with distance. Since too large an in- 
crease in 3 with increasing distance is unrealistic and since Y "1 implies 
that the temperature is constant, we used an Intermediate value for Y of 
1.05. The hydrodynamic solution then gives a value for the Initial radial 
velocity at the solar surface of 8.2 km s”^. The initial meridional 
velocity is assumed to be zero. 

The magnetic field is assumed to be force-Iroc initially. The 
equations for the vanishing of the Lorentz force and the absence of free 
magnetic poles then determine the initial magnetic field configuration, 
which is given in terms of Legendre polynomials of order 1. The initial 
configuration selected is the dipole configuration obtained by using the 
Legendre polynomial of degree 1; i.e., 

3 B cos 0 

B » , 

r 3 ’ 

r 

B sin 0 

■> 0 - 

r 


where B^ is the reference magnetic field. This configuration tends to 
restrain material ejected near the equator and consequently is referred to 
as a "closed field" configuration. The reference magnetic field is 

selected to obtain the desired value for 3. Values for the reference 
magnetic field of 0.83 G and 2.35 G at 1 at the equator yield 
values for the plasma beta of 4 and 0.5, respectively. 3 increases with 
radius and decreases away from the equator. The initial value of 3 at 
1 at the pole is a factor of 2 less than at 1 R^ at the equator. 
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Th© above chosen Initial state does not, of course, represent a 
steady-state soCution to th® complete two-dimensional equatiisns. Hence, 
when used as the initial state, the time-dependent solution to Rquations , 
(1) will evolve until the solution relaxes to or approaches n steady- 
state solution (the coronal streamer). 

The remaining quantities that must be specified for the numerical 
solution are the boundary conditions at the four boundaries in Figure I, The 
solution is assumed to be symmetric about the equator, and hence, the 
boundary conditions there become (using the notation in Figure 1) : 



P ^ P« 
1 3 





For the coronal-streamer calculation, the pole can also be considered as 
a symmetry boundary. The boundary conditions there then become the same 
as for the equator with the following exceptions! 





Note that there is no flow across either of the above two boundaries. 
There is flow, however, across the remaining two boundaries at 1 and 
5 R (the inner and outer boundaries) and this makes their treatment 
somewhat more difficult. The reason can best be explained using the theory 
of characteristics. This has been done previously for one-dimensional, 
hydrodynamic flow by Nakagawa and Steinolfson (1976) and Stelnolfson 
and Nakagawa (1976). We will briefly discuss the extension of their 
analysis to the preseiit problem (see also Endler, 1971). Equations (1) 
are a hyperbolic system of equations, and therefore, the characteristics 
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for tlie equations are all real. In thia eaRu Inrommtlon from any tnitut; 
la the flow fitUcI can only propagate In the eharneteriatle. aireetUnvs. In 
the radial direction, then, information propagate.a at tlve ain veloelUieB 
equal Co Che sum of the radial velocity and the radial component of each 
of the six chnracCcriscic vcloclclea. At the outer boundary the flow la 
ganorally supersonic and super-Alfvenic, and hence, information from the 
outer boundary only propagates downstream, i.e., outside the region of 
interest. Another way of viewing this is tl\nt all .six radial charaoteristle 
directions at the outer boundary are pusltive. This lmplle.s that the 
boundary conditions at 5 It can be specif led arbitrarilv» We ehoo. 4 e to 

III ‘ 

linearly extrapolate vising the relation + AQ (notation as in Figure 1) 

where Q in any dependent varinble. and Aq is calculated at each time .step from 
the values at the two radial grid points adjacent to the boundary. At the 
inner boundary two of the six radial cluwactorlstlc directionB are negative, 

and consequently, information from the region of interest propagates up- 
stream to the boundary. In this ca.se four dependent variables at the lower 
boundary can be specified arbitrarily, and two miust bo calaulatevi from some 
fornrof compntibiUCy relations (see Stelnolfsou and Nakagawa, 197b). 

Strictly speaking, the compatibility relations are equations that can he 
derived from Equations (1) wlvich must be, satisfied by the dependent variables 
in each of the cbaracceriatlc directions, Steinolfson and Nnkagawii (1976) 
have shown that first-order or second-order (linear) extrapolation often 
works as well as using the more complex compatibility relations. We choosy 
to hold the pressure, density and radial magnetic field at choir Inltia)/;' 
values at the boundary and calculate the meridional velocity so that the 
total vc.locity and magnetic, field are parallel at the boundary. The 
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radial velocity and meridional magnetic field are calculated from linear 
extrapolation. When the velocity approaches ^ero inside the closed-field 
region, the pressure and density at the boundary are allowed to increase 

V\ 

accordingly by linear extrapolation. 

B. Numerical Results 

The initial-boundary value problem discussed above forms a well- 
posed problem, and the time-dependent equations are now Solved numerically 
using the appropriate initial state and boundary conditions. It remains 
to select a criterion by which to determine when a steady-state, coronal- 
streamer configuration has been achieved. Admittedly this decision must 
be somewhat subjective. The criterion we used was to continue the time- 
dependent solution until the configuration did not change appreciably 
over a period of at least one hour of physical time. In general, the 
smaller the value of 3, the more rapidly tlie solution evolved to a coronal- 
streamer configuration. 

1. 3 “ 4 

The initial dipole magnetic field is shown in Figure 3(a). 

The field lines are spaced at 5 deg Intervals at the solar surface. 

The curve consisting of long dashes represents the location where the 
radial flow velocity is equal to the sound speed (sonic curve) . The 
radial flow velocity is equal to the Alfven velocity based on the total 
magnetic field along the curve composed of short dashes (Alfvei\ curve) . 

The evolving magnetic field is shown for three successive times in Figure 3. 
The numerical solution is changing so slowly with time after 24 hours 
that the solution at that time is selected as the final coronal-streamer 
configuration. Several of the dependent variables after 24 hours are 
shown in Figure 4. The nondimensional values given for the pressure and 
density are the values referenced to their respective initial values at 
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each point in the flow field. From Figures 3 and 4, it is easy to see the 
field lines evolving from a closed dipole field to a coronal-streamer 
configuration with the closed field lines lying beneath and adjacent to 
open field lines. The sonic curve is displaced inward in the final state, 
except for a small region around the equator, due to the general increase 
in velocity, The velocity decreases somewhat at the equator fro"' that 
in the initial state and increases at the pole as shown in Figure 5. The 
changes in the pressure and density at the equator and pole referenced to their 
respective initial values at each point are also shown in the figure. 

The top of the closed-field region in Figure 4(a) at the equator occurs 
at approximately the same radius at which the Alfven curve intersects the 
equator. The reason for this is that the magnetic field at the equator 
is approximately zero above the closed region, and the velocity is approxi- 
mately zero in the closed region. There is one field line above the inter- 
section of the Alfven curve with the equator that appears to be closed in 
Figure 4(a), but this is believed to be due to the interpolation between 
the relatively large grid spacing that we have used. The Intersection 
of this field line with the equator is gradually moving outward. 

The pressure and density are increased over their initial values 
and the velocity is approximately zero in the closed region as indicated 
in Figures 4 and 5. The temperature is proportional to the ratio of the 
pressure to the density, and hence, it is evident that the temperature also 
increases throughout the flow field. The increments in pressure and 
density between the contour lines in Figure 4 are equal so more tightly 
grouped curves denote larger gradients. The velocity vectors closely 
parallel the magnetic field lines in Figure 4 despite the relatively large 
value of 3. The open field lines never become completely radial J not 
even the field lines next to the pole. Hence, assuming that Che flow 
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I 

i, ‘ft 

does follow the field lines, the flow near the pole is in a region where ^ 

the croBS-sactional area of an infinitesimal flow tube increases more 
rapidly than as the radius squared. Several steady-state studies for flow 
in coronal holes where the cross-sectional area of a flow tube is specified 

T 

to increase more rapidly than as the radius squared have appeared in the 1 

literature (see, e.g., Kopp and Holzer, 1976, and Steinolfson and Tandberg- 




Hanssen, 1977). These steady-state studies conclude that the velocity 
increases in coronal holes, but that the temperature and density decrease, 
which is in contrast with the present results. The reason for the 
discrepancy is not immediately apparent, but it is most likely related to 
the different boundary conditions applied at the solar surface in the 
two studies and to the liiclusion of the meridional terms in the equations used 
in our study. 

The results discussed above for 0 « 4 display several of the character- 
istics expected of coronal streamers (correct magnetic field configuration, 

higher pressure and density and ze^o velocity in the closed region, etc.) 

"\ 

in spite of the fact that 3 = 4 is 3)arger than would be expected in coronal 
streamers. We now consider a reference plasma beta almost an order;of 
magnitude smaller (3 = 0.5) which is a \^alue more consistent with observations* 


r 

I 

I 

I 

I 

i 

i 
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2. 3 = 0.5 

The Initial dipole magnetic field configuration and the location 
of the sonic curve are the same regardless of the value of the reference 
magnetic field as can be seen by comparing Figures 6(a) for 3 = 0.5 and 


I 

I 

I 


3(a) for 3=4. The only difference in the two figures is that the Alfven 
velocity Increases with the magnetic field causing the Alfven curve' to 
move outward as 3 is decreased (magnetic field Increased). As the magnetic 
field evolves in this case, as shown in Figure 6, the Alfven curve moves 
inward throughout the flow field as does the sonic curve except for a 


I 

9 
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small region around the equator. .‘>ome o£ the variables in the final 
coronal-streamer configuration are presented in Figure 7. The sonic 
curve innersects the equator in Figure 7(a) at approximately the same 
radius as for $«4 in Figure 4(a) indicating that the sonic Mach number near 
the equator must be similar for the two cases. The radial velocities 
along the equator are about the same as can be seen by comparing Figures 
5 and 8, and therefore i the temperatures along the equator must be 
similar. The temperature distribution along the pole is clearly quite 
different for the two values of 0 with the temperature for B-0.5 
being considerably larger. The velocity in the open region near the pole 
for 3 - 0.5 is more than a factor of two larger than the velocity over the 
closed region near the equator. The velocity in the closed region is 
again approximately zero as would be expected, but the pressure and 
density in the closed region are increased considerably more for 3*0.5 
than they were for the larger value of 0. 

The distribution of 0 throughout the flow field is shown in 

Figure 9 for both the initial state and the coronal streamer. 

For the initial state 0 increases from 0.5 (0.25) at the solar surface 

to 10.2 (2.6) at 5 R along the equator (pole). In the coronal streamer 

0 

the value of 0 Increases at the equator and becomes infinite along the 
neutral line above the closed region. Since 0 is so large near the 
equator and since the incremental spacing between the curves in Figure 9(b) 
is fixed, the 0 distribution is not shown within 10 deg of the 
equator in order to show the distribution in the remainder of the flow field, 
0 Is decreased at the pole near the solar surface and increased at 5 at 
the pole compared to the corresponding values in the initial state. Through- 
out most of the flow field, 0 is lower in the coronal streamer than for 
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Che initial state. Because of this and since Che thermal pressure is 
higher In the coronal streamer than initially, the magnetic pressure in 
the coronal streamer is also higher than the initial value throughout 
most of the flow field. 

As for the previous example, the velocity flows more or less along 
the field lines. Also, none of the field lines becomes exactly radial. 

The radial location of the top of the closed field region is not changed 
appreciably from that for the higher value of 0. The dimensions of the 
closed-field region are considered below for a wider range of 3, 

3. Coronal-Streamer Structure 

As the magnitude of the reference magnetic field is increased 
in studies such as the above two, one would expect that the field lines 
would remain closed to a higher altitude in the solar atmosphere. That this 
indeed happens is illustrated in Figure 10 where the dimensions of the height 
and the base at the solar surface for the closed region are plotted as a 
function of 0. The circles indicate values obtained from individual computer 
runs. A closed region is not achieved for 0 *= 100. The base dimension 
of the closed region is relatively unaffected by the values of 0 unless 0 
becomes larger than about 5. The maximum values of the pressure and 
density (referenced to their respective initial values at each point) in 
the closed region are plotted as a function of 0 in Figure 11. The large 
increase in the pressure and density with decreasing 0 is due partly 
to the increased Lorentz force which is directed inward on the closed 
region and partly to the increased magnetic field strength which Increases 
the tensile strength of the field lines and makes it more difficult to open 
the closed field lines* Tiiis latter effect has been demonstrated by Steinolfson 
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Qt al. (1978). Both of these effects must be balanced by an increase in 

thermal pressure in the closed region, and both effects also tend to 

compress the fluid in the closed region. For 3 * 100 the magnetic field 
is so weak (fl - 0.17 G) that the solution remains essentially unchanged. 

The velocity and Che thermodynamic variables at 5 R are presented 

© 

as a function of meridional angle in Figures 12 and 13, respectively, for 
several values of $. The velocity is consistently larger over the pole 
than the equator with the velocity differential between the two regions 
Increasing dramatically as 3 decreases. The pressure and density are 
highest near the equator and the temperature is highest at the pole for 
each value of 3 - The temperature distribution for 3 equal to 4,10 and 
100 are so close to being the same that the curves cannot be resolved in 
Figure 13. For the two lowest values of 3 (0.5 and 0.1), the density 
decreases and the velocity increases from their respective initial 
values in the polar region in analogy with flow in coronal holes (Kopp 
and Holzer, 1976). However, according to the study of Steinolfson and 
Tandberg-'Hanssen (1977) the temperature and pressure should also decrease 
in coronal holes, while in the present study both variables increase. Hence, 
although the open magnetic field in the present study is diverging as 
would be expected in coronal holes, it appears that the flows in the open 
regions in the present study are not analogous to previous steady-state 
studies of flows in coronal holes. As mentioned earlier, the discrepancy 
may be due to the different boundary conditions in the two studies, but 
it may also be partially due to the meridional gradients in the equations 
that we use which are neglected in the steady-state studies. 
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IV. CORONAL TRANSIENTS 


The coronal-streamer configurations discussed In the previous section 
are long-lasting structures which may remain virtually unchanged on the 
solar surface for periods up to several days. These configurations are often 
observed to exist prior, both spatially and temporally, to the occurence of 
mass ejections or, as they are often referred to, coronal transients. We now 
simulate a coronal transient by using the final state for the coronal streamer 
as the Initial state for the coronal transient. 

A. Initial State, Boundary Conditions and Perturbation 

As mentioned above, the final state for the coronal streamer is 
used as the initial state for the coronal transient. The reference 0 
referred to in this section (and all other reference conditions) is the 
reference value used for the coronal streamer. Naturally, all the dependent 
variables may have a different value in the final streamer simulation than 
they did in the corresponding initial state. 

The boundary conditions applied at the equator, pole and 5 R for the 

0 

coronal transient are the same as those used for the coronal streamer, At 

the solar surface a perturbation in the thermodynamic variables from their 

initial values (the final values for the coronal streamer) is introduced 

near the equator to simulate the solar event that produces the resulting 

coronal transient. The perturbation is assumed to occur instantaneously, 

to be symmetric about the equator, and to extend in the meridional direction 

5 

over a total of 10 deg (5 deg on each side of the equator) or 1.22 x 10 km. 
The perturbed values are maintained for the duration of the calculation. 

The magnitude of the perturbation varies with each simulation and is discussed 
along with the corresponding numerical results. The pressure and density 
axe assumed to be constant throughout the perturbed region. Outside this 
region they are calculated from linear extrapolation as for the streamer 
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simulatloiVf The remaining variables at the solar surface are simply set 
equal to the value at the grid point next to the boundary at each time 
step; i.e., Qjj ■ (first-order extrapolation). 

For the streamer simulation the computation had to be 
carried out for the entire flow regime at each time step. In this 
section we only compute the solution from the solar surface to the 
leading edge of the coronal transient, 

B, Numerical Results 

The coronal transients for the two values of 6 used In the 
previous section ($ • 0.5, 4) are Investigated. The transient for the 
lower value of 3 is examined in more detail since, as will be shown later, 
it agrees more closely with observations. Some of the more important results 
for other values of 3 are also discussed, 

1. 3 ** 4 

The thermal pressure in the perturbation is increased by a 
factor of 2.6 above the Initial value at the solar surface for the streamer 
calculation which, as discussed in the previous section (see Figure 5), 
remained unchanged throughout that calculation. The pressure increase 
results from a factor of 2 increase in temperature and a factor of 1.3 
increase in density. Some of the variables in the resulting coronal 
transient are shown in Figure 14 after 60 minutes and in Figure 15 after 
120 minutes. The corresponding variables for the initial state for the 
transient are those in Figure 4. The time is now referenced to the start 
of the transient simulation. The pressure and density are referenced 
to theLr initial values for the streamer calculation at each point in 


the flow field. 


Th« lavge velocity that developa at the aolat aucface In the vicinity 
of the pectucbed region causea the flow to become super-Alfvenlc there 
after 60 minutes and both auperaonlc and auper-Alfvenlc after 120 minutes 
as seen In Figures 14(a) and 15(a). The coronal transient Is preceded by 
an MHD shock which has the characteristics of a fast shock since It 
increases the meridional component In both the velocity and magnetic field. 

It Is not clear from the figures that the shock normal (or the direction of 
shock propagation) lies directly along the roognetlc field, except, of 
course, when the shock propagates Into the open region along the equator. 

Note that a meridional velocity component cannot be produced by the shock 
along the equator because of the boundary conditions there. When the shock 
normal and magnetic field are parallel for some finite length of the shock, 
the fast shock becomes a swltch-on shock since it produces a meridional 
component in the magnetic field and velocll'r when neither was present ahead 
of the shook. The numerical shock is necessarily spread over a few grid 
points due to the relatively large grid spacing that we used. Hence the shock 
location shown in the figures Is selected to be approximately at the mid- 
point of the jump in the dependent variables produced by the shock. The 
shock has not completely formed at large meridional angles near the pole 
after 60 iidnutes. The shock velocity along the equator increases from 306 km s 

at 2 R to 584 km s”*" at 4 R Implying that the shock travels considerably 

0 0 

faster in the initially open region than It does in the Initially closed region 
One reason for this Is that there la an ambient outward flow In the onen 
region which tends to increase the shock velocity in the laboratory frame. 

The contact surface separates the initial coronal plasma from that which 
has been emitted from \’he perturbed region. The maximum increase In pressure 
lies along the equator below the contact surface while the maximum Increase 


in dnntlty U b«tw«tn th« contact nur£ac« and tha ahock. Due to the relatively 

larga valut of 8 umd for thla almulatlon, the inagnatlc field doea not have a 
large effect on guiding the outward-propagating plaama and la, to a large 

extent) eaaentlally carried along with the plaama. Aa can be aeon In the 
flgurea* the contact aurface travela outward both radially and meridionally. 

It appaaea from the flgurea that the fluid la propagating acroaa the field 
llnea (which ahoiiild not occur due to the Infinite conductivity) alnce the 
contact aurface at 120 mlnutea encloaea more field lines than at the earlier 
time. That thla la not actually happening la a result of the method used 
to calculate the field lines. The field llnea are calculated by starting 
at the same location on the solar surface at each time step and using the 
magnitude and direction of the magnetic field at that time as f|he Bcarting 
values. Since the magnitudes of both components of the magnetic field at 
the solar aurface are allowed to change with time, different field lines 
may be traced out from each starting location at each time. 

2. /9- 0.5 

The thermal pressure In the perturbation Is Increased by a factor 
of 10) the temperature by a factor of 5, and the density by a factor of 
2. All of these Increases are with respect to the corresponding final 
values at the solar surface for the streamer calculation. Some of the 
variables In the coronal transient produced by this perturbation are 
shown In Figure 16 after 80 minutes and in Figure 17 after 180 minutes. 

The corresponding variables for the Initial state for this transient 
are given in Figure 7. 

The relatively low velocities and large magnetic field prevent 
the flow from becoming either supersonic or super-Alfvenic at the solar 
surface near the perturbed region as It did for ^ ” t. The transient 


i» igAlit pr«c©dttd by a fast MIID ahock which may, ovar a poetlen of its 
OKCary;, b« a swltch-on ihock. The ehock velocity along the equator 
increaaes from 397 km at 2 R»to 517 km at 4 R*. Hence, ai for 

the reiulti at «« 4, the ahock tr»/ei8 fatter in the initially open 
region than in the initially cloeed region. However, the ahock velocity 
is now larger at 2 R* and emaller at 4 R§ than it waa for /9» 4. The 
ahock velocltlea for these two valuea of 0 and for other values of 0 
will be discussed In more detail later, 

The contact surface does not travel as far in the meridional direction 
as for 4. This is due to the larger magnetic field in this case 
which tends to restrict the meridional movement of the plasma. Along 
the equator the maximum Increase in pressure again occurs below the 

contact surface and the maximum Increase in density between the contact 
surface and the shock. This is shown better in Figure IR which gives 
the radial profiles of the radial velocity, pressure and density along 
the equator Initially and at two subsequent times. The mass excess 
initially In the closed region of the coronal streamer (illustrated 
by the increase In density between 1 and 2 R# in the figure) is at the 
later times distributed between the shock and the contact surface. 

The fluid ahead of the contact surface (Including the fluid Initially in the 
closed region) has been heated, compressed and accelerated by the shock. 

This figure also demonstrates how the shock is spread over a number of 
grid points. At t - 80 minutes, the shock, which should be a sharp 
discontinuity, extends approximately from 3.75 R« to 4.25 Rq. 

The curves consisting of the dark shorter dashes in Figures 16 (d) 
and 17 (d) enclose a region in which the density is higher than in the 
adjacent regions in the meridional direction. That is, in a scan from 
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the equator to the pole at a fixed radius, an Increase in density would 
be found in this region, this is illustrated by the meridional scans 
of pressure and density at 4 Bn for three times in Figure 19. The 
density has a definite peak at approximately 30 degrees at 180 minutes 
where the scan is now through the "leg" of the high density region in 
Figure 17 (d). At the previous times the scans are partially through 
the top of the high density region. The pressure, on the other hand, 
does not have a peak away from the equator at 180 minutes and decreases 
more or less monotonically for all three times from the equator to the 
pdle. The curves with the dark longer dashes in Figures 16 (d) and 
17 (d) represent the approximate locations of the maximum in the hlgh- 

denslty region. The Importance of this dense region will be discussed 
in the following sections. 

The distribution of /9 throughout the flow field (except near the 
equator, since /9 again becomes very large there) for the two times in 
Figures 16 and 1 7 are given in Figure 20. The Initial distribution 
of /3 is shown in Figure 9 (b). The /8 contours are very similar to 
those of the pressure — especially after 180 minutes when the magnetic 
field is almost radial. These plots Illustrate that for the region 
shown, /3 remains relatively small during and following the passage of 
the coronal transient. 

3. Other Pertinent Results 

The shock velocities along the equator at various values of j3 are 
shown in Figure 21. A wide range of perturbations in the thermodynamic 
variables was used to obtain these results as shown in Table I. The 
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«hoek vtlocitles foi‘ the 4 and 0.5 atudiea diacusaed abova aca 
r®pi*e8an(:ed by tha dashed curves. Theiro are aaveral ljuportant reaulta 
that should be noted from the figure. Pirat of all» despite more than 
an order of magnitude change In both /3 and the magnitude of the pressure 
perturbation, the shock velocities for all simulations do not vary that 
much In magnitude} at h R^, the velocities are all within 200 km s * of 
each other. Note also that at a fixed a reduction In the magnitude of 
the pressure perturbation reduces the shock velocity although the shape 
of the curve remains generally the same. This is true for both /3“* 4, 
where both runs have the same density perturbations, and /9» 0.5, where 
the two runs have different density perturbations. On® would expect 
that If tlm density perturbation was varied over a wide range this 
result may be altered. Finally, as decreased, the shock velocity 
curves tend to become flatter until for |3 « 0.1 the shock velocity Is 
essentially constant over the radius range considered. Tf the results had 
been extended to a larger radius, the curves for 0*4 would certainly 
approach a constant velocity also although the final velocity In this 
case may bo somowhat larger. 

The fact that a larger pressure perturbation was used to obtain the 
results In Figure 21 as 3 decreases is not coincidental. In fact, at 
each value of (3 there is a value of the magnitude of the pressure per- 
turbatlon (at a fixed value of the density perturbation) below which a 
physically realistic solution (In the sense that tlie flow along the etiuator 
is outward) is not obtained. This would be expected since the perturbation 
is in a closed-field region, and a certain pressure force is required to 
open the field lines vso the fluid can flow outward. As the strength of 
the magnetic field Increases (13 decreases), one would expect that a larger 
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forde would bo roqulirod. The approximate minimum prosaura perturbation 
for a fixed density perturbation of 1*3 is showi\ in Figure 22 for the 
range of /9 used for the results in Figure 21. The circles indicate the 
runs used to obtain the curve. Hence, it is rlear tliat as la decreased 
below 0.5, the required pressure perturbation increases very rapidly, 
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V. COMPARISON will! OTHER APPLICABLE MODELS fj 

The other models which we will consider here are those which are similar j 
to the one we have developed In the sense that they consider the transient to 

i I 

be the coronal response to a simulated solar event near the solar surface. 

That Is, the driving force responsible for the transient originated and remains ^ , 

i 

near the solar surface. We do not consider the magnetically-driven model 
proposed by Mouschovlas and Poland (1978) and further studied by Anzer (1978). 

They assume that the driving force Is the magnetic force produced the 

' i 

Interaction of currents flowing In the loops with their self-induced magnetic 
fields. Yeh and Dryer (1980) have attempted to show that coronal loop transients 
> cannot be driven By such magnetic forces. 

All similar models that have appeared in the literature which purport 
to be applicable to the formation and propagation of coronal transients 
in the lower corona (with the Inner boundary inside the critical points 
near the solar surface) in either the meridional or latitudinal plane 
and including the magnetic field have assumed that the atmosphere is 
initially stationary. On the surface this appears to be a reasonable ‘ 

assumption since the kinetic energy in the inner corona is certainly less than 

. V 

either the thermal or magnetic energies. It is definitely simpler to assume 
an initially stationary atmosphere since an analytic planar solution Including , 

i ■ 

the magnetic field (other than a purely radial field, which Is unrealistic In 
the lower corona) and 'j./ finite velocity has not been obtained for the lower 
corona to the author's knowledge. In addition to including the velocity in 

I 

the initial streamer for the transient in the present study, the magnetic field 
in the streamer decreases with radius similar to a dipole field. This offers the | 


advantage of producing a reasonable Increase in 3 with radius. For any hlgner 
order field 3 Increases too rapidly with radius to be reasonable. We first 
present some results for an initially stationary atmosphere and then compare these 
and other similar results from the literature with our results obtained 
using a coronal streamer as the initial state. 

A. A Model with an Initially Stationary Atmosphere 

The initial state for this model is taken to be an Initially hydro- 
static atmosphere for the thermodynamic variables and a dipole magnetic 
field. This initial state represents a solution to the steady-state 
form of the complete set of time-dependent equations (Equations (I)") . The 
reference conditions for the thermodynamic variables and the magnetic 
field are the same as used for the streamer simulations. Again we will 

consider the /3 = 4 and /3 “ 0.5 simulations. The boundary conditions, 
grid spacings and spatial distribution of the perturbation are also 
the same as used above. Hence the simulations considered here are in 
every way as close as possible to the same as those done for the coronal 
transients in the previous section. 

The spatial distribution of some of the variables for the coronal 
transient for ^=4 are shown in Figure 23 after 120 minutes. The mag- 
nitudes of the pressure and density perturbations used for this simulation 
are the same as those used in the previous section for P- 4 . SitBilar 
results for /3= 0.5 are presented in Figure 24 after 80 minutes. The 
pressure in the perturbation is now increased by a factor of 6 and the 
density by 1.3, which is not the same as those used in the previous 
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section for /9 » 0.5. The corresponding results obtained when the initial 
atmosphere Is a corona] streamer are given In Figure 15 for ^ ■ A and In 
Figure 16 for 0 *■ 0.5. The pressure and density along a merldlrnal scan 
at 4 Tv0 for the present simulation with 0 • 0,5 are Illustrated In Figure 
25. These results are analagous to those In Figure 19 for the streamer 
Initial configuration. All the above results for the initially stationary 

I 

atmosphere will now be discussed vls-a-vls the results from the previous 
section. 

B. Comparisons 

For convenience, the results obtained in this section and other 
previously published results for an initially stationary atmosphere will 
be referred to as the "stationary” results, while those in the previous 
section for an Initial atmosphere consisting of a coronal streamer will 
be referred to as the "streamer" results. The figure numbers of the results 
presented in this paper will only be referred to when some confusion may 
arise. 

Since the perturbations are Identical for both the stationary and 
streamer results for 3 “ A and since the magnetic field in this case is 
relatively small, a comparison of the results should be similar to what 
one would expect from a similar hydrodynamic study. Indeed they are, as 
will soon become apparent. The shock travels more rapidly and the fluid 
velocity is larger for the streamer results due to the initial fluid 
velocity. The shape of the shock curves in Figure 23(a) and l5(a) are 
also similar, which would be expected since the shock velocity in a lab- 


oratory frame depends to a large extent on the characteristic speeds in the 

ambient medium through which it is propagating (in addition to Its dependence 

on the velocity in the ambient medium). As the magnetic field is relatively 

small in this case the sound speed dominates, and, of course, this speed is 

approximately the same for the two sets of results. The separation between 

the shock and contact surface is larger in the stationary results since the 

velocity of the contact surface depends entirely on the fluid velocity which 

is smaller for the stationary results. The general shape of the pressure 

and density profiles for the two cases are quite different with the higher 

pressure and density values being confined to a region closer to the equator 

for the streamer results. The maxima of the pressure and density are larger 

for the streamer results because of the Increased values in the initial 

»* 

closed magnetic field region near the equator and the solar surface. 

The magnetic field would be expected to have a larger effect on 
the results for /3 ■ 0.5. Note that a larger pressure pulse was used for 
the streamer results when /3 * 0.5 than for the stationary results. This 
difference is not expected to affect the comparisons. The shock travels faster 
along the equator (despite a smaller pressure perturbation) and slower at the 
pole for the stationary results as compared to the streamer results (Figures 
24 (a) and 16 (a)). This result is reasonable since shocks with a given initial 
strength always travel faster across field lines than along them. This effect 
is not important for /3 = 4 since the field is so small. The other results 
mentioned above for j9 = 4 remain much the same at the lower /9. It is 
Important to note that the high pressure and density region near the equator 
for the streamer results is even more exaggerated at /3 = 0.5. 


The curves conslstlrtg of heavy short dashes in Figure 24 (d) enclose 
the high density region that would be observed in a meridional scan at a fixed 
radius In analogy with the slmillar curves In Figures 16 (d) and 17 (d), A 
curve with heavy longer dashes again Identifies the approximate center of the 
maxima of this hlgiv-denslty region. The density peaks away from the equator 
for the stationary results show up on the meridional scans at 4 In 
Figure 25; however, they are much further from the equator than for the 
streamer results in Figure 19. 

The maxima of the high density region (the curves conslf^tlng of the 

heavy longer dashes in Flgureslb (d), 17 (d), and 24 (d) are shown for 

several times for /3 ■ 0.5 and for both the stationary (dashed curves) and 

streamer results In Figure 26. The maxima for the streamer results are 

contained within 30 deg of the equator, while the maxima for the stationary 

results extends to the pole. A plot slmillar to that in Figure 26 is shown 

in Figure 27 for only the streamer results for j9 ■ 0.1. This plot demonstrates 

that as P is decreased below 0.5 the angle of confinement of the density 

maxima decreases. The variation of the density along a meridional scan at 

4 R for 0.1 corresponding to the results in Figure 27 is presented in 
® / 

Figure 28. This result is in comparison with that in Figure 19 and shows 
that as P decreases the density peak becomes more prominent. The importance 
of these results will be discussed in comparing the results with observations 
In the following section. 
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VI. COMPARISON WITH OBSERVATIONS 


Among the several classes into which coronal transients in general 
tend to be grouped (Hildner^ 1977), the type we have attempted to simulate 
is the loop-like transients in which mass is ejected from near the solar 
surface to the outer corona. The major observed characteristics of this 
type of transient will first be discussed followed by a comparison with 
the calculated results. We concentrate on attempting to simulate just 
the characteristics for this single class of transient since it is 
extremely difficult to simulate all observed aspects of one particular 
transient (due to Inhomogenletles, non-planar effects, dissipative effects, 
etc.). One would also expect that different physical processes are 
Involved in the different classes of transients. 

A. Observed Characteristics of Coronal Transients 

The velocity of the leading edge of a transient is a relatively simple 
quantity to determine from the white-light observations. These velocities 
generally lie in the range of 100-1200 km s ^ with an average value of 
470 km s”^ (Gosling et al. , 1976). The velocity tends to increase or remain 
about constant between approximately 2.5 - 5 R . 

The reason that coronal transients are observed by the white-light 
coronagraph is that the density in the transients exceeds that in the 
surrounding medium. Hence it is necessary that the simulated transient 
produce a density enhancement similar to that observed. It should be 
pointed out that the coronagraph only observes the projection of the 
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fcranslent: in the meridional plane eo It la certainly posaible that they 
do not lie In that plane. For the purpose of the following discussion 
we will assume that the observed transients do lie In or near the meridional 
plane. As shown by Hlldner (1977) > the observed transient tends to take 
the shape of a loop, with the loop extending from near the solar surface, 
out to some maximum height (which Increases with time) In the atmosphere, 
and back to the solar surface. The northern and southern edges (the legs) of 
th^ loop generally remain stationary so the loop just stretches outward 
with time. The centers of the loops usually lie within 40 dog of the equator 
with a symmetrical distribution about a maximum at the equator. The lat- 
itudinal extents of the transients are usually less than 65 deg. 

The coronal magnetic field cannot be measured directly so other 
means are necessary to infer Its value. Dulk et al. (1976) used data 
from type IV radio bursts along with white-light data to calculate the 
magnetic field for the loop-shaped mass ejection of 15 September 1973. 

They calculated a value of 2.6 G at a height of 3.1 R« and a 3 at that 
height of 0.007. Further evidence for expecting relatively large fields 
near transients Is due to the observations that at least three-quarters 
of the transients occur over active regions (Hlldner, 1977) where the 
magnetic field Is larger than elsewhere on the sun. Mass ejections also 
tend to be associated with magnetic phenomena (flares, eruptive prominences, 
etc.) which occur near lines on the solar surface which separate uni-polar j 
regions of the photospherlc magnetic field. 

The mass associated with a mass ejection generally lies in the range J 
1 - 24 X 10^^ g with an average of 6.2 x lO^^g while the energy ranges ” 

from 1.6 - 69.8 x 10^^ erg with an average of 11.9 x 10^^ erg (Hlldner, 1977) 
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In order to obtain these values it is necessary to make some rather 
speculative assumptions concerning the size o£ the transient, 

B. Comparisons 

The radial velocity of the leading edge of the transient can be 
easily simulated by the model proposed herein, as shown In Figure 21, 
or by any of the previous models which do not use the streamer conflg'> 
uratlon as the Initial state. Larger velocities than those in Figure 21 
can be obtained in our present model by Increasing the magnitude of the pressure 
perturbation. The observed increase In velocity with radius, however. 

Is generally between that shown In Figure 21 for - 0.5 and /9 - 0.1, 

Implying that the lower 3 results are more realistic In this sense. 

The fact that previous models also reproduce the observed velocities 
(as demonstrated by the one-dlmenslonal study by Stelnolfson and Naka- 
gawa, 1977, and the two-dimensional study by Dryer et al., 1979) makes 
the leading-edge velocity a poor quantity by which to judge the relative 
merits of various models. 

One observed characteristic that previous models cannot adequately 
simulate (at least not for any of the results published to date) is the 
shape of the loop as given by the density enhancement. This is seen 
quite clearly for the calculations presented In Section V for a non- 
streamer, stationary Initial state (see Figures23(d) and 24(d)). The 
enhanced density region moves both radially and laterally away from the 
simulated solar event similar to the motion of an expanding bubble. Similiar 
results for a non-streamer, 'Stationary initial state have been obtained in all 
previous two-dimensional studies (see eg. Dryer et al. , 1979; and Nakagawa, 
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Wu And Han, 1978). Thaa* pravlotii rtauXt* art obtainad ragardlaas of 
vhathar tha Inlfclal magnatlc flald configuration la aaaunad to ba cloaed 
over tha location of tha aimulatad aolar avant (aa foe tha initial dipole 
configuration for tha atraamar calculation in Figure 3(a)) or open (aa 
if tha configuration in Figure 3(a) wan rotated 90 dag so the configuration 
originally at tha pole waa at the equator). The edges of the loop do 
not move laterally in the present model as seen in Figures 16(d) and 
17(d). Tha moat graphic comparison of the present model with previous 
models is shown in Figure 26 which shows the outward motion of the high 

density regions near the outer edge of the transient. The latitudinal 

I! 

extent of the transient for the present model is 60 deg which lies within 
the observed latitudinal extent of 65 deg. The results in Figure 27 show 
that if 8 is decreased the latitudinal extent of the transient decreases. 
Hence it is obvious that the present model agrees considerably better 
than the previous models with the observed density enhancement which 
Is certainly one of the more important observed characteristics. 

Since the coronal magnetic field cannot be measured directly, a 
direct comparison cannot be made between observed and calculated magnetic 
field magnitudes at various points in the flow field. The Inferred value 
of the magnetic field of 2.6c at 3.1 Rn made by Dulk et al. (1976) 
agrees closely with the values at that height in our model for the 8 » 0.1 
calculation. Dulk, however, obtained a 8 of 0.007 at that height which 

is about an order of magnitude less than what we calculated indicating 

:s|l' 

a difference in thermal pressure. Although the magnetic field cannot 
be measured, it is reasonable to expect that if the 8 ia relatively small 
(as observations indicate it should be) it may be important in guiding 
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th« lluld; that it, In deterinlnlng th« shape of Che expanding loop. 

Thus It is logical to assume that the magnetic field prevents the lateral 
expansion of the loop, as has been done by several authors. This assumption 
would agree with our model since It Is the magnetic field that contains 
the lateral expansion In Figures 26 and 27. The fact that mass ejectlona 
generally originate near lines separating unl>polar regions on the solar 
surface Implies that previous studies using open field lines do not apply 
to maaa ejections. 

We have not attempted Co calculate the maos or energy associated 
with the mass ejections In our simulations. Due to the assumptions that 
must be made for the dimensions of the transient both in the observations 
and In the simulation, it Is not felt Chat a direct comparison of the 
absolute magnitudes would be meaningful. 


VII. WSCU8SI0N AND CONCLUSIONS 


The miiln objactlve for d«v«Ioplng self-consist;<int models such ee the one 

» 

we have prederifced herein la to attempt to gain a better underatandlng of the 
phyeical mechanisms occurring in and responsible for observed coronal 
transients. As such, it is necessary that the simulations reproduce the 
available observations as well as possible in order that valid conclusions 
regarding the behavior of transient quantities Chat are not observed can 
be made from the simulations. From the previous section it is apparent 
Chat the present model does reproduce the observations better than similar 
models published previo\>sly. Again we should point out that the only ocher 
models we are considering here assume that a sudden disturbance near the 
solar surface produced the transient. 

The major difference between the present models and previous ones 
is the Initial state for the coronal transient. While previous models 
simply use a stationary atmosphere with a force-free magnetic field, we 
have used a coronal- streamer configuration in which the atmosphere is 
not stationary (except in the closed-field region) nor is the magnetic 
field force-free. Coronal streamers are frequently observed in the lower 
solar corona and are observed to occur over uni-polar regions of the 
photospheric magnetic field. As mentioned previously, this is the region 
in which mass ejections commonly originate. We should also mention that 
the class of transleifits which we are concerned with is the loop-like 
transient in which mass is ejected to the outer corona. Although it does 
not appear likely that previous models apply to these transients, they 
may certainly be applicable to other classes of transients or other solar 
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phenomena t 

The magnetically-driven models studied by Mauschovlas and Poland 
(1978) and Anzer (1978) may also merit some consideration as the mechanism 
Involved In mass ejections. However, the simplifications made In their 

I 

analysis and the fact that their study Is confined to just the top of 
the loop limit the application of their results. A more complete study 
of magnetically-driven models may prove worthwhile. ' 

tn the present and previous similar models, the radially- 
propagating, leading edge of the transient Is a shock wave. This 
shock wave certainly agrees with observations In that It travels radially 
outward and produces the necessary density enhancement. The fact that 
both type II and type IV radio bursts (which are Indicative of shocks) 

are often observed In transients substantiates the poss:|Lbility that 

// 

the leading edge may be a shock. Other than the frequent observation 
of radio bursts there is no a priori reason why the leading edge should 
be a shock. In fact, in the magnetically-driven models. It is not. 

The lateral edges of previous models are' also shock waves (which travel 
laterally as well as somewhat radially), while In the present model 
they are due to plasma confinement by the magnetic field. Since the 
lateral edges of observed transients do not travel laterally, fehey, 
of course, cannot be shocks. ii 

The available observations imply that the magnetic field associated 
with mass ejection coronal transients would be expected to be relatively 
large. The plasma beta is certainly less than one and may be one or 
two orders of magnitude smaller. This conclusion being due to the usual 
occurance of mass ejections over active regions and their association with 


magnotlc phenomena. The present study supports this view since the best 
agreement with the well-observed, leading-edge velocities of transients 
is obtained for ^ S 0.5. The lower beta results in the simulation 
also produced better lateral confinement of the edges of the transient, 

As P Is decreased, the more the simulated coronal transients with 
a coronal streamer initial state differ from those with an initially 
stationary, force-free initial state - as shown in the results in Section 
V, One reason for this Is that - as discussed In Section III - the 
coronal streamer diverges more from the initially stationary, force- 
free state with decreasing /3 . For the coronal streamerj as ft decreases 
the following changes occur in the variables: (1) the pressure and density 

In the closed-field region Increase, (2) the velocities Increase throughout 
the open region and the gradient In the velocities between the values at 
4 R (5 at the equator and at the pole Increases, and (3) larger changes 
occur in all the thermodynamic variables and, as for the velocity, the 
gradients In all the thermodynamic variables between the values at 4 
at the equator and at tlm pole Increase. Hence the merldloi\al inhomo- 
geneities In the streamer increase with decreasing 

For a constant reference temperature, as we have assumed, the only 
other parameter In our model besides /3 is the polytropic index Y* We 
used a value of 1.05 which Is net much different than the constant temp- 
erature value of 1. From the results in SectioitsIIl and IV it appears 
that the temperature does not change much in our simulation of either 
the streamer or transient providing /3 > 1. Hence for a large /3, the 
constant temperature solution should not be much different than ours. 
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However , as /9 is decreased below 1 , the temperature changes more and 
more from the steady-state value Indicating that a constant temperature 
solution is no longer valid. 

The mass that is contained in the expanding loop in our model (i.e. ; Che 
enhancement that would be observed by a white-light coronagraph) is entirely 
composed of mass which was in the corona prior to the occurence of the solar 
event responsible for the transient. Part of the mass was initially in the 
closed region and part in the open region of the coronal streamer. In this 
respect our model contrasts with other more approximate one“dimensional 
models ichere the mass is assumed to initially reside in a loop without any 
consideration given to how the mass-filled loop was formed (Mouschovias and 
Poland, 1978; and Anzer, 1978). The transient in these simplified models consists 
of an expansion of the top section of the initial loop outward through the 
corona, while in our model both the formation and propagation of the loop 
are calculated 'self-consistently throughout the meridional plane. 
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FIGURE CAPTIONS 


A schematic of the portion of the meridional plane In which the solution 
Is calculated! The grid spaclngs are also sketched (not to scale) 
along with the nomenclature used for the boundary conditions. 

A schematic of the magnetic field lines In a coronal streamer. 

The evolution of the coronal magnetic field for /9 "4. The sonic 
curve Is represented by the curve with long dashes and the Alfven 
curve by the curve with shorter dashes. 

Planar maps of several of the dependent variables In the coronal streamer 
for ^ ■ 4 . The maximum values given are the maximum vaJ ues throughout 
the regions shown. The velocity vectors point in the direction of the 
velocity at their base and their length la proportional to the magnitude 
of the velocity. The pressure and density are referenced to their respec- 
tive Initial values at each point in the flow field. 

Radial distributions of the density, pressure, and radial velocity 
initially (solid curve) and in the coronal streamer for /9 » 4 (Figure 4) 
at the pole (shorter dashes) and equator (long dashes). 

The evolution of the coronal magnetic field for /9 ■ 0.5. See the 
caption for Figure 3. 

Planar maps of several of the dependent variables in the coronal streamer 
for ^ ■ 0.5. See the caption for Figure 4. 

Radial distributions of the density, pressure, and radial velocity 
initially (solid curve) and In the coronal streamer for 0.5 
(Figure 7) at the pole (short dashes) and at the equator (long dashes) . 

Planar maps of the plasma beta for the /9 «« 0.5 simulation initially 
and In the coronal streamer. 

The dependence of the dimensions of the closed region on the reference 
beta. A closed region is not obtained for /8 “ 100. 

The ildependence of the maximum pressure and density in the closed 
reg:jion on the reference l^eta. 

Thcj' meridional distribution of the radial velocity at 5 R as a 
functlop; of the reference beta. 

The meridional distribution of the thermodynamic variables at 5 R as a 
function of the reference beta. ® 

Planar maps of several dependent variables in the coronal transient for 
j9 ■ 4 after 60 min. The approximate locations of the shock (double 
lines) and contact surface (heavy dashed lines) are indicated. See 
the caption for Figure 4. 


15. Planar maps of several dependent variables in the coronal transient 
for /9 <■ 4 after 120 ntln» See the caption for Figure 4. 

16. Planar maps of several dependent variables In the coronal transient 
for 0.5 after 80 min. The short dashed curves In (d) enclose 
the maximum density region and the curve with longer dashes traces 
out the approximate location of the maximum density In this region. 

See the caption for Figure 4. 

17. Planar maps of several dependent variables in the coronal transient 
for /3 <■ 0.5 after 180 min. See the captions for Figures 4 ‘and 16. 

18. Radial distributions of the density, pressure, and radial velocity 
initially and at two later times in the coronal transient for/3" 0.5. 

19. Meridional distribution at 4 of the pressure and density in the 
coronal transient for 0.5 at three times. 

20. Planar maps of the plasma beta for the coronal transient for^* 0.5 
after 80 and 180 min. 

21. The variation of the shock velocity along the equator with radius 
for several values of the reference beta. 

22. The magnitude of the pressure perturbation required to produce a 
mass ejection coronal transient for a fixed value of the density 
perturbation as a function of the reference beta. 

23. Planar maps of several dependent variables in the coronal transient 
for /9 “ 4 after 120 min for an Initially stationary atmosphere 
with a force-free magnetic field. See the caption for Figure 4. 

24. Planar maps of several dependent variables In the coronal transient 
for 0.5 after 80 min, frr an Initially stationary atmosphere 

with a force-free magnetic field. See the captionsfor Figure 4 and 16. 

25. Meridional distributions at 4 R, of the density and pressure at three - 
times in the coronal transient for /3 ■ 0.5 for an initially stationary 
atmosphere with a force-free magnetic field. 

26. The maxima of rhe high density regions at several times in the coronal tran- 
sient with /3=“ 0.5 for an initial state consisting of a coronal streamer (solid 
curves) and an initially stationary, force-free atmosphere (dashed curves). 

27. The maxima of the high density regions at serveral times In the coronal 
transient with /8» 0.1 for an Initial state consisting of a coronal 
streamer. 

28. Meridional distributions at 4 R of the density at .two times In the 
coronal transient for /9 =• 0.5 for an Initial state consisting of a 
coronal streamer. 
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APPENDIX A. COMPUTER CODE FOR NUMERICAL SOLUTION 

The simulations of the coronal streamer and the coronal transient 

■1 

are performed In two separate computer runs. The numerical solutions 
for the streamer simulation at each value of 3 are output to mass 
storage devices at selected time Intervals. These solutions can then 
be recalled by the plotting routine In Appendix B and plotted In 
various formats by the computer. The streamer simulation at the 
last time for which It Is calculated Is then used as the Initial 
state for the transient simulation. The transient solutions are 
also output to mass storage devices at selected time Intervals, and 
plots can again be produced. 

As an example of the computer time required, the streamer simu- 
lation for 3 ■ 0.5 required 6 min and 45 sec on the CRAY-1 computer 
to simulate 16 hrs of physical time. The transient simulation for 
3 * 0.5 required 40 sec on the same computer to simulate 2 hrs of 
physical time. 

Since the codes for the streamer and transient simulations are 
essentially Identical, only the code for the streamer simulation 
is presented here. The variables that must be input to the code 
are given first followed by a listing of the code. 
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r“ I** 

5a 

Variable ColuEms 


Description of Input Variables for Computer Code tfflPHS* (CRAY-1) 


Definition 


Coffiaents 


Units 


m 

Title 


R0 




1-4 

5-80 


1-8 


Run Identification (Card 1) 


Run number 
Run title 


Right-adjusted integer 
Any alphanumeric information 


Steady-State Parameters (Card 3) 


Reference radius 


a. The time- dependent solution is computed for 

r 2l R0 

b. The values of the remaining variables on this 
card must be the values at this radius 


solar 

radii 


T 

9-16 

Temperature 



5 

D 

17-24 

Electron (or proton) density 



no/cm^ 

PI 

25-32 

Polytropic index 

Constant everj-where 
set to 5/3. 

(1< PI< 5/3). If input as 0, 


B 

33-40 

Magnetic field 

The reference value 

at R0 and theta = 90°. 



Pulse Parameters (Card 4) 


TAU 

1-8 

Pulse duration 

DMG(1>, 

9-16 

Determine 

variables 

R0. 

the dependent 
in the pulse at 

I = 1,9 

17-24 




* 

73-80 




minutes 


The values i-, .ut depend on the type of pulse being 
considered - aee IPULS on the option card 


This code is two-dimensional and is applicable to the meridional or r-r plane in spherical coordinates. The azimuthal 
components of the velocity and magnetic field are neglected. This code does not use the large core memory. 






Variable Columns 


Comments 


Units 


REND 

DR 

1-8 

9-16 

DTT 

17-24 

DRT 

25-32 

TMAX 

33-40 

P^MAX 

41-48 

DTTP 

49-56 

DRTP 

57-64 

TMAX 

65-72 

DTP 

73-80 

DIPT 

1-8 

STBY 

9-16 


Definition 


Numerical Solution Parameters (Cards 5 and 6) 


- Card 5 — 
Steady State 


Maximum radius 

Grid spacing in radial 
direction 


1 solar 
radii 


Grid spacing in theta 
direction 

Grid spacing in radial 
direction 

Maximum distance in theta 
direction 


Time Dependent 

It is best to select DTT such that R*DTT is about 
equal to DRT at (SMAX+R0)/2 


The solution is counted from the equatorial plane 
(theta - 90“) to TMAX >90*. 


Maximum distance in radial 
direction 

Theca increment at which 
solution printed out 

Radial increment at which 
solution printed out 

Tine limit 

Time increment at which 
solution printed out 


The solution is computed from R0 to RMAX. 


If either (or both) is input as 0, every grid point 
is printed out in that (or both) direction at each 
specified time increment (DTP), 


The solution is calculated imtil the simulated tine 
exceeds 

If input as 0, the solution is printed out at every 
time increment. 


— Card 6 - I 

Tine increment at which If input as 0, the solution is written out at every I 
solution is written out on a time increment J 
plot file 


Stability constant 


a. STBY :< 1.0 

b, RecoiSBend using a value as near 1-0 as possible 
without numerical instabilities developing. 


degrees 

solar 

radii 

degrees 


solar 

radii 

degrees 

solar 

radii 


minutes 



Defiaitioa 
Smoothing constant 


Variable Coluans 
SMTH 17-24 

TSTR 25-32 Times at which solution 

written on mass storage 


t 



« 43 


Ccinr-its 




. Removes some numerical instabilities. 

. Recommended using SHTH 1.0 if possible- 
c. Only used if ISHTH « 1 and ItSCK « 0. 

The solution, can then be used as the initial state Hitaites 

for restart in a separate prograss. 


! 



Variable Coliuans 


Definition 


Value 


Description 


Ootions (Card 2) 


ISTST 

1-4 

Steady-state solution 

1 

IPDLS 

5-8 

Pulse 

2 


1 


2 


ICHCK 9-12 Tine-dependent solution check 0 

1 


Boundary treatment at R0 


Isothemal solar wind. 

Poly tropic solar wind. 

The pulse depends on the particular study being 
conducted. The types of pulses currently available 
are the following (all sysmetrical about theta = 

90 "): 

Pulse in tenperature and density. 

I9tG(L)=TM (Meridional extent of pulse in degrees). 
DMG(2)*6T, MlG(3)*5p. The naxiraua pulse in T and 
p is at theta * 90* , and it falls off sinusoidally 
to ambient values at theta = TM. See IPDLS = 2 also. 

Pulse in the aagnetic field. DMG(l) as for IPLTS*1, 

DMG(2)*(Bg) gaTiss, DMG(3)*(By.) where (Bq) and 
P * P P 

(By) are the maximum values of the latitudinal and 
^ P 

radial components, respectively, in the pulse. Input 
0 if no pulse desired in one of the cc«^onents. The 
theta dependence is the same as for IPDLS=1. For 
both IPULS=1 and =2, M^(4)=DTI where DTI is the 
time in minutes for the pulse to reach maximm 
magnitude. The pulse increases linearly in magnitude 
with time. 

Do not check the time-dependent solution. 

Check the time-dependent solution by starting it 
with the steady-state solution. The space terms are 
handled properly (and the steady-state solution is 
correct) if the time-dependent solution does not 
deviate appreciably from the steady-state solution 
after a sufficiently long time. 

The boundary treatment depends on the type of pulse 
(IPDI^) being used. See the subroutine BDSY. The 
current options aret 


IBDRY 


13-16 



Stsady-state magnetic field 


calculation li'i 



IPTFL 


IPTVR 


ISTSL 


ISMTH 


XTMSX 




IFVCSTMTrs. 21-24 


25-28 


29-32 


33-36 


37-40 


41-44 


45-48 


Contact surface calculation 

* 1. 

Print out radii at which sound 
and Alfven Mach nuabers are 
initially and at each tiae 
print out solution ^ 

Dependent variable plots 


Store solution at tiae TSTR 
on mass storage for restart 

Smoothing 


Smoothing 


Boundary treatment at 
90*. 

KP^iS 

:v 


No boundary treatment 


V and Vg at R8 set to value at adjacent radial grid 


point (R#<-DRT). 

^Closed dipole 
0 

Open quadmpole 
P Open dipole 
Closed hexapole 
Open hexapole 
'Closed decapole 
Open decapole 

If IMAGF - 1, 4 or 6 (2, 3. 5 or 7), the option 
IBDX0 is set to 2(3) intemallv. 

b 

Follow the contact surfaces along r at 8 « 90* 

■h a ^ 

Do not follow the contact surfaces along r at 6*90* 

Print out jH 

No print 

No plots of the depe^eot variables 
Produce plots of the dependent variables (UAH) is 

X! I) ™ ■ 

Produce plots of the dependent variables (VCAR) “ 

fci « 

Do not store solution* ^ ^ ' 




“r*” 

No saoothing applied. ■ 'vl'lW . 

Lapldus saoothing aethod applied. | 

The input value sinus uae is the number of tine cycles 
that are skipped before sroothlng. Only used If tSMTH- 
X • 

The nosenclature is: < value a: ^0-DTT, t ^2 " 



value at 90+DTT. For all cases 


“3’ ’^l * "''3’ 


*ei L*93’ ®rl 


=•3. Pj 


C=! 


nrz>r'-"t 











IBDXM 


IBDYM 


IRDTN 


49-52 


Bounr!ary treatment at TMAX 


4 


0 

1 

2 

3 

53-56 Boundary treatment at RMAX 


57-60 


0 

1 

2 

Radiation-Not included in 0 

these codes 

1 

2 


IL01JB 


61-64 


Application of numerical 
differencing procedure 


0 

1 


Mo boundary treatment 

The nomenclature is: ( “ value at TMAX-DTT, ( )2 * 

value at TMAX, t * value at TMAX+DTT. 

No boundary treatment, calculation stops when reach 
TMAX. 

Qj = ^ 2 * where Q is any dependent variable 

= 2 Q 2 - where Q is any dependent variable 

where Q is any dependent variable 

The nomenclature is: ( * "value at RMAX-DRT, ^ 2 * 

value at RMAX, C * value at RMAX+DRT. 

No boundary treatment, calculation stops when reach 
RMAX. 

= 20^ ~ Q dependent variable 

where Q is any dependent variable 
Neglect radiation 

Calculate radiation (Cox-Tucker), but don't allow 
it to effect the solution 

Calculate radiation (Cox-Tucker) and include it in 
equations. 

Use usual method with equations in conservation form. 

Use modified method. This method is often necessary 
w’hen 3 is small (3 ~ 0.1) and the magnetic field is 
essentially radial. 
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JOB NO. MHDHS Input Format - 2 5 April 


NO. 

2 3 4 5 


:R;u>i I 


;i.V^CL^ j 

. R.N;...' 


xT- 


.f. 


COST CC*~RZ. NC 

Fortran STATEMENT * 


lABa 


t 

’*4 


7 3 5 !0 n n 13 34 !S J4 17 1» 17 20 21 22 23 24 25 26 27 28 77 3C: 3» 32 32 34 3S St 37 33 34 40 41 41 4J 44 i4S 4* 47 *S 47 iC J1 32 S3 54 S3 34 57 5* 61 *2 63 44 44 *6 t’ 4* *« 79 7' 72] 73 34 •• 7* 77 7* 7» tO! 


dUjn »ti f ip a t lanj ;(|Ca r dj Ij) :j . , : . | ^ 


T~ 


it it l!e ' 


t • . - ■ - -|— 


OE^tl o|i 4 S_,- All right-adjusted integers (Car d 2): t - 

|iSTSt] IljULs| 1 CHCk[ IBDRy| IMAGf]i5WCS j IPTFl] IPTVr[ ISTSL | ISMTh[ NTlilSM| IBDXo] IBD3i^[ IBDYM | IRDTO ( n.0^ 


S t en dy 1 _ jS t a t e 

Pa rame te r s 
• pV- ' ■ 

( Ca rd . 3 ■) : 
1 1 


1 ■ 

■* * 

i ' R 0I • ‘ f T 

, . i. . ■ 


1 PI L 

B 

1 



FtI 'Is e ! IP ar ;am e t er s (fc a rd ; 


T AUi 


% 

] Dt^Ts) I BHG(6)^ D!^(7) J »»G(8) [ DHG(9> j 


! 


:iu_me ^r]i ;ca 1 , So^lu t in n' Pjax^ me t e^rsj j; :Ck_rd , 

RRHD’ I [ DR I _ : ; D R" ; | i 

--/yr, \T. ^ 


-4 4 ' ♦ * 


miuX f DTTP 1 DRIP i T«AX 


I TMiX 




^ DTP ^ 


N ua eir 

L ‘ D 

4— ..-4 


ic ai Solu'tio'n jPa r am et eyjs ?(;Cja rii 6.).:^ ; j . j 




ST, BY ^ ^ BMiT'H' 




X?S?T R^ 



I . 4 . ^ ^ - 4 . f- — ♦ „i' - 

9 -» *-- * -•4'-' "I -i ■ '4 '■■ ■•' *■ 

, . ' * i ■ 




: 
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Listing of Cot ■*uter Codes MHDHS 


S7tglc Qi-rZS/iO 05’.-3St19 C->iJ 

ro *jo'‘,si 43 ,;*A 7 i-'ci,s-piCri stei%,’^f$on, mhohs 
rr3TTT?^ LlTr6‘*=f^«i TcT,'?,.- 

THIS CCjE is ;f*PL:C^oL" TO TPc ’'E9IDI0MRL PLA«E IN SPHeSICAL 


icira»:fajMCfai 


>■»! 




(X,Y5 =LA^E. THE TEP»tINOLQGY li - X FEF=RS TO TKHTA CT) , T TO R, 

O ~G VT, V TO VP., '■^TO ?T, ‘ EE. 1X 7 TO T, IJT? TO_ 

ct:*--' jv xcixTj,y(iY7J ,'Uix7,i>',Z) ,v(ix7 , iy7,2).mix’', iy7,2J. 
1PCI»7, IY7.2),lD 1S'-Y (1X7) ,3Y(IX7,IY7.7J, SXf 1X7.1 Y7,2),A25tlX7,IY7 


2AIZ 

COrPOV /COaMI/ YO,TO,00, TC,AO,RttO,PO,S&,TIO,YENO,DXT,DTT,X«AX, 
1Y«AX,TAU,T.XAX,DTF,KS,V0,l.'O,':O.ArO,AREAS, Af eat 

/(.V'-'w:/ L ST L Tx , » 7 TtMticT rPGuT;'>t«p ,»>«p,'ffTT,'STBYTP 1 7i> p P‘ 

CC -'-ON /Ce-ILX/ IS:-:TP,5FTB,.-.~' I.'.ILOW.; 

CI'--'"CY /-".••rd[>/I"DXO,IrSXA>-,IE:.Y'' 


CCP»3V /COhPOY/lAEOY 
DIPHASIC?. TITHT) 


2 FOC-^iT U<.,TV»i / 1 ?iA) 

3 FCF.-';i (• iOPUT -• // 3X*RUN IS rVTI FI CATICW ' / SX’rUW AUKF EP 


4 FCn‘ 


SX'EO =*,1PE1Z.S,' SOLAR RAOII’/SX'T =’, 

, • <'/EX'0 =',E13.5,' rLECTSONS/CHj* f Sx*"! =*,DFn2.'/5X 

: TTUYX'S 5 ^5 

r r 7 Ei-.v) 

T { / 'X'PLLtS car4-£t 5-~' / aX’TAU = ’ ,T FE 1 0. 5 , ' *)I-.UTES*7 
= ',En. 3 ,"' t.tC'FecS* / 3 > 'OT =', Ell. T/rr^ro ■ = •■■; cTrrj/Tx 
=',Eld.3,' “IHUTES*) 

T (ITE-^.O /3FH..'.) 

T~r7~Tx'^S’"l7&Y-5f AT?' * TjOf.'P ARAPETjcRS'" / 5x* Rfc" L =' - iF 
• SCLA" =4111' / 3X'D- =’,E13.5,* SCLXF RA'<li* // :-x 
-LEPELOiJ.T ECLPTIV.’ =>-'i ETEPS' / EX’OTT ='.E12.E. 












no 

"0 

CO 


R£n> ( I \ , S 3 Tor'' * \ ,'t i t T T c ^ ^ “ o l s , i c « ck T i" *> c r yTi*" * >? R V i f *• c s' 

1, IRTFL.IOTV:, IRTVv,:S->n,‘.-T - sy, I.UXC , :3r>”, I=?C.Y'-,ir-TK 
?, IL '-. v . 3 , T-fOY _ _ - 


If IRTV- = 2 


s'Ao (v-,A> 

oRirt (i;,;3Yc,7c,30,Fi,:-c 

2 F Ck2 •£^«0«03 ^ I — - »/ 5« 

REA3 <IV,t> TAU,X?'P,W''P,C''R.V»-XP'>,Y»'XP2,TAU‘> 
CO TGCIY, 12,‘5),1FULS 
^ ^ i 1 ^ ^ T3 » F *j f X'* ^ ) V* P ^ • V ^ X P J 

2Z PFiM n, re-, ,y^t ,H- p,: -? ,Y xPF 

W*’p = «-'!r/rO 


o» f=dr;p/ oc 
GO 


i.:, xm'Vp 'P.'2"p, ’ 


10 RFA" Cl',, ) V£PP ,cr , ;. /T,l Y *, X -AX, V^S.< fDXTP.rY'P ,J"AX,; T_P,D JPTJL,^ 

ISTGtl'STo'i' -- - 

AP.ITC C3C, / lYpf.i ,I Y ,GXT,&YT,X':ftX,YKAX ,vX7P , ”YTP ,T vs y , j.-TR ,0 YRTL , 
ISTpy, 3 

eO'Lt^TF'OT^' t IE’St , : Y .'J"ACf , » '£ TOStPi .iPXf» ,3 PT-.P,.. VTP } 

C/Ll T ;”v.R - C Tv' i. Y, '»-:p,OY''r-, Y:-'*p’1 ,W Xp ' ,"S ■■'I ,3 R,,; . rSApXE , ARE AC, 

"'i‘TF£, ,v fOC,*. ^ ^ . 

cfTo? " " ',_j— 

E'l:- 

$i : F;‘.rTM= r:^'orp c 3t-.sv,£#.’p ,>.- yte ,ykxpi ,Y'’Y 07,-i-.«’ ,i r«c£, area" , 

■»•-''■•••;-.---«■-• . p ,- ^. -. - --,-*^.2 ■•■•£. e "* '"'' 

rC'TCCi. T;’v'-pprC%ri\T iOttiCC' CalCULATIO*: 

= Mfc-t.-ETP'. < :x-’='-, iy- = s?: 

— "cT'«ow'xTn'TT,Y{TY?>TMi23rY,iY'^,;i r.’c ;x?; jv"».T>, ''irx’^.iY'-,? 3, 

1PCIX7, :-,'?,3 5,Li:c"Yc:x?j,.»YC:y.?,:Y7,2),''x(:x7,r''i',C)«Ac' (I x’’,iY?', , 
:A7Z CIX?,IV7) 


E : ■ F N i : ; V : a a 2 1 : y . i y - j 

11 cr * i ” HT '/3 4' '7 p '***^*'*, **«*,♦, ,,»»*,*****,-* / *" x '^*'* TI 

3 o * 7 3Cj, 3'^P********"**'*'*’'****’'*******^ * ^ 

Tr*7tf;»? f^Tr<-ri ar- = « , 


'- CiR’cHaSAiT 


:;» 'T E'T 


yCCx 


, ; , p t . , 


t X r X • - 


Ty Y“' j x'.;' 7 ‘ " 3'i, > ; ‘ 1 * * • ? * X 

16 rC-"'.''T U 7 , :t C 1 P -: 1 / 

'3 3’-:''T £;<’«p\px:v 


■3 I ^ t -* ' IPA-I r -,' ^, i.p »,«'-• Y , 




4 i ■ 


o o to o o p o o L> o o p o o 
oopoopoopoopoo 


II 


. •'•.AL =’,E1'!.5T'. fClXETIC 5«»V IT»T10XAL =-‘,'fcii"5 

e/ 5X»MASS 3i)lL0-UP =V,E11.5,’ {.PS/kI^Z =*,E11.5t* PAtTI CLES /<’'! * > 

? FOpi'AT tSX'F CnaSACTEKISTlC OlPECTiONS AT RO* / 1 AX *T • 1 OX 'LJ' 

1 1 0X ■ LA ‘ 1 u X ■'l 5'^f’' * ' 1-6 *■> ~ ” 

2,t)j 

-RO =‘,1PE10.A," sOlAR RAtll =*,E10.A, 


YC0NTS=Y (1 5 
YCTSR= YC1 ) 
LCOM-1 
LCOr.TI =? 
LCT?=1 


N-PTVRsO 

IF (IPTVft.EO.I) ’JOTVR = 1 


PLLoO=0 
TF = XCc)-1 
TP = 


5RA=6. 67E-0 
SUfJF=1 .yf ?E?2 


0X02=0X1/ 


tXTT2=CXT-»H 


TTK,lL«X,IFFXrIFPY 


otftl=otptl/tio 

TI"PL=DTPTL 


CALL LftX 
1' (iPCFY.HE.O) 

1 CALL 30RY tlP.KTIvc ,lJ6iiY,lLVX,Ti»',IFP.L'*XPy) 

'call OrOpY 1 1 Lr ? ri =r,P , 1 >pVi P5"’ ~ ' 

IF (ICHCK.EO.P) r TO 1 
1L'-X = 1 


00 












W t> 


2 LpISt>Y(l) = ' 

C^LL H?53Y(IP) 

i 

I Jail, f CLiifT" <: IT^ t , ”\ :T, y'-ViiVyvxF ^ y 

If (IP7’Y«. £ -.0) .-O TT : 

J[F£LSC=1 _____ _ 

IPULit=? 

liFITE uc,n) 

-FiTi * ■■■ ■■ ■" ■ 

^7:? ; = ?,Ti‘ e + 1 

I = C:cHC‘-.t-.1) •- TP 2; 

1 T C 1 £ N D 1 . T)~cr lC~ FuT'?' i Af : £“ 

1Y'-XP2, L'‘X = r ,L''XCY1 ,TA!Jn 
Z: IF (Tl-.LT.T'iAX) SO TO ' 

£ 1= (ICf.C< .E '. T3 C-T TO 1 
IF C.TI-J.tT. 1) 
iCAi.!. iC!SYi UP, 

:r CM £,.<T.i5 CALI or-TiYiTiMY.iEM 
“it CALL LAKY ( 1 = ^ . It , £ 7T«C , XL“X . I FF X, I FF Y T , I f PY2 ) 

z'F V:oh:k ov'" .1 > ’call Fs TYi c iF> 

>.TT»/0. 

OC, T t 0 = T ,L'' Y - 

DP t y I = StL'X-'1 ■ 

IK-I-1 

o’ <i'-, j)+oTr» T»{vf j, : p3««x t :, o .i 

" E ■' c 1 ,7', I p > I , j , I - p 5 * M c i ; j , : ► F 3 -0 c i';o , : F p > 

FCOLTO CPOTACT ivOFACiT 'tONO T AT X=^ 

3 F COFyOO -.-S3 00 it 

Y c e •• T s = Y c c T ; * 1 T T y * c V '< r, L c O'- T , I ff >~+ c V c; •-• T 3 - -rc 1 
1IFf 3-V <£ ,LC”*.T,rHP}3/(YTLrOM13-Y(t_COT,T3 3> 

_IF CYCCtTr.LT.YfLrOfiTlJ) b'T TOJi? 

L'iMI =L:r\Y-^T 
'? :f cti • . 1 . 0 ,Y„f5 r- rr. i>. 

Y C TSn= YC TS f TottmIITCT; LCIpTTpf YC fs'nM cf^ 

1 VCS ,L€ TA , iPS')) / { Y(i.t IP-! 3 -y ltC'!«3>) 

IF CY£TA.y.i.T.YC;.crS13> ;.C TO it 
LF 

L£Tp-3 = LCTF + 1 
ff::.t OtT -^ALOFi 

? IF ( T '"FTamTiI'N M 
Tir" = Tir.»Ti" 

wOITc <10,123 TIK.Tr.‘'',MI->,F,0T7. 

ir"fl!F»tr.T- .'U~iT‘‘!Z I: ■ 

YClM'-sYCpO’S-Yf 

'< C T =YCC*. TP» Ty - / ^ c 5 

Y C'T S «= YCTTTAYn ~ ~ 

y:t3p>;=y£t:f»p3/i ,H5 


XCY-S ,TAU1 3 


~c Y ,TI“" , IP ,T>* XPl'," 


Y C ',y ,I yP> ; 


crt.T3T*Fr, ti" lccSt 


3 *ivi ‘,L£ TPt ,IPP3 








f\i r\i 


fcBLT=0 .0 
E3K 1=0.0 
£SL<=0 .0 


D3U1=0.0 
DeUZ=0 .0 



P^!S1=3 (I.lrIPPJ 
CO 22 J=2,L-"Y 

b I . -J=i^ .'r ( 1 , J , I P- J ^VTI i J ,Tr7*'^ 2 

T“=:==>n,j,if-) 

Sr£c=DCI,J,IP}*Ci;Cr,J,lP)**2-»VCI,J,I») 


F."S3=D (It J, IP) 

Ef s'2=Eij?«a*0YO2*(bT21+:T?3) 
EM 2=t fr J ^ ■♦■"0’Vw ? 1 T 1 VTH t 

esk;=elki+dyoz» cf kei+pkeZ) 
ErGZ=£EG'*OyOZ* Ct»!r1*GREZ) 


0 

0 

0eL2 = DLLc*CY02* CSrS1 + PK32) 
TH£1=THE2 

il 

\t t.'c. 

*■. LI = y 
£ 1 =o r ? 2 
«r*:si = ftws2 

%) 

e 

<• 

£PLh=EtLJ'*DXO<;» 1 EE '.) + £■:>=* 27 
EFLT=EELT+0XC2*CE.-T1 + EcT 2) 
E£.LSf = £r L»:*3ye2-<E*" K l + £aK2) 
















3j. K j ce= citj C8- c r i j . 

32 RJC3)(=-ST(I,J.IPJ*RJCB 
RJCeY=&XCI ,J, rP)*RJC5 


Tr TjTbTTTj uO 7t» 

IF (ICHCK.es. 1) GO TO 1° 
TEKPR = PCltJ t IPP>/0CI, JflPP) 


• 611(1. J 


' KRiTt (10.14) J.yrro.ud.J.iFP) tVCi.J .iPP).BJcax, dci.j.ipp). 

, ip-<i.j,ipp) . axfi.j ,iPPi.eyti,j,iPP).«.icBy.TEHPR. bet* 


bO la cu 

19 TEKPR=P(I.J ,I‘'>/D(I,J .IP) 

B£TA=eETAO»PCI, J,IP)/(6X(I, j.IF)*»2*eY(I,J.IP)»*2^ 


' [p ^5 >#4 , Bx ( I . j »a.p ) ,Bx 1 1 7 j ,V> j Ci WeiPP #pctA v 

"'” 20 YIP =NTP*eVTP 


uitfmiti 




iriUE 

IF (I.GT.c.OR.IRDTN.EC.C) GO TO dO 
FF.IOT 62 


. DO 61 J=1 ,L£UP 

>. IF CYCJ).UT.rTP.AN0.J.N£.l.SU9) GO TO 61 


T ii’ K= h' 1 1 , j . j, r J # u i , j . 

QK!)M = R ADFBCTE"PR»TG)»(0(2, J,1P> *00) «• 2*CHA0N‘ 

aKD.NM=CRD*J-rAOTK'(J) 

■>'=55f^,ei;<’.eONTIN,U6'- ■■ . 

r. ^60 yTP * rCl)-1.E-6 •' 


XTP = XTP+OXTP 

17 continue 

IF CIPTFL.cC.1) go to 86 


PRINT 95 


RSE1=Y0 

KA£1=YG 


IAE1=1 

RFS2=va,1,L)/S'PTCPCI,l.L)/oa,1.L)) 


Rr-t2=V(I,1 ,L)*S 0 r.TC 0 CI ,1,L)/(eX(1.1,0**2-‘BY(I, 1,L)**2)>/A20 
IF CRMi2.CE.1.) IAE1=Z 


60 TO C8V,9f>) 
S9 SFS1=SFS2 


IF CR-1S2.lt. 1.) SO TO 9C 

ISH1=2 


90 GO TO (91,V2),I*c1 
9^ ESA1=RyA2 

— T^^?7?r=TrrT7rrT3 E r c 0^ I riTE^Tf? X Cl , jTO rt c* r Y ti . J , L V • 

IF CK-'’ A 2 .L 1 .1 .) CC TO S 5 


DBMmmssami 


yiaaw 


92 I? (lSt1*£S-2.A''D.IA£1.E(..2) GO TO 93 
3S CONTI'JLE 


T'T 








v.n V/t fS I' ^ ^ 4> ^'l^w W v4 l/i W W Ul Ul ol M M M <M M W f>J r\> *4 -k ■».* 

\^rv>J ^ c> << ou -Ni cH ui -r* OJ w <a H o* vn w w i-o ->| o *o itl -si *> /4 4^ W fv -* C3 >« cd *>* 




93 XKM=XC D/BAPDG 

PRINT 9«,I,XKH ,RSE1,RAET 

=7 COKTISlOE _ 

.'6 XTF=x(c.)-1,r-6 
TlfP = TI'"P + OTP 

CRfftTr PLOT fILtS 

40 IF (IFTVH.EC.C) go ' TO Si 

IF (TIP.tT.TIRPL.A?iD.Ich3*Ea.1J GO TO It- 

56 TIKMP=TIJ«*TIO 

LClSRl =tula>X^ 

CO 57 J=2,L!)1S^X 

57 tPI YFI CJ)=L0IS,NY(J)*1 

60 To (6o,6A>, IPTVR ' ' 

66 ENCODE (12,i7,lH02) TIHPP 

WRITE OO) <lHD2{I),I = 1,I)t UDlS*1,<t.DlTB'J(J},J=2,t01SnX) 

CTT'To ci "■ ■" ' ' ' 

ft. WRITE (10) TIRf'P, LOlSW-1 , (lDI YP1 1 J) ,J=2,L9ISKX) 

A5 DO ia 1=2,LDISRY 

L0IS«2=U0IS’'Y(1)»1 
GO TO (66. £>7) , IPTVR 

66 WRITE (10) CP(I,J,IP),D(I,j,IP),SX(l,J,IP3,aY(l.J.IP)»UCI,J,IPT, 

Tv ClTr, IP ) , "■ ■jS1,L01S'-2)' 

GO TO 4? 

67 APC1 = X(I) 

IF (Il*AGf .EO.T.oF.IXAGF.EC.S.Oii.lwXGfZta^'^ AfcGl=SftG5^ilC5 

00 55 J=1,LDISf<2 

59 A7Z(I-1,J> = A2F(I-1,J)*Y<4)*E1*»(ARG1) _ . _ , 

"lifiTtT~rrT5~rp ( I7 j",i r ) ro <i , j » i p> ,=x < tvo . i pj . ^Vti , j , I p> , u < I . J , IP) . 

IV (I, J,IP),A2ZCI-1,J) ,J=1,L0ISK2) 

ts ccrtinoe 

MPTVii^HPTWRfi ~ ; 

PRINT A5.NPTVR 

IF (IPULS6.£t.1) GO To 56 7 ... 

TT:' PL= Tl .-.RL+CTTrL ■ 

51 IF (I£ND.Ei.2) RETURh 

to TO 16 - 

eiTD : ^ ^ ^ 

FORT PAS, 0==>LIS,BM=STSTC1 .Ft 

SUfrROOTINE STST USTST , DY , IPASF , ASEAg ^AR EAD ,IPT Ft., IPTVR , OX TP> 
cal colat c stgABT -TTATr^TOrrioN- T^;^iiL''Tii>-r-!>EFST<ingNr'ARBrr5 ' 
PABAnETEB (1X7=39,1X7=51) 

CCW'CtJ Xdx'^) rY(lV7).U(IX7,lY7,2),V(IX7,:Y7,2>,C(:xT,;T7,2), 

TFnx77TY7,2),LDis«Taxt'i',smx7-,is7,2j,pms7,rYr7?J7i2?TrxT7rY?> , 

2AZZ(;IX7,IY7) 

CC)»*«0N /COM**)/ YOtTO.DO, lO.AO.RSO , PO, SD.TIO, YESB ,» XI ,6Yt,X«*X» 

TYmTTw;t»TxvBTF;5?Tvo,«;F-o,A^Apfc«^ 

CCR'’0N /CC«’‘27 L0ISWX,L''X,L*'Y,ICHCK,IFBLS,XPP,.WP,OTT ,ST«Y ,PI,J**P 

COK(’?OK /C0»RDK/mDTH,RA0T?>(IY7)»CRAB«,«ADTVt<Ix7,IY7) 


1 FOF’‘AT ( // ' OUTPUT -• // 1CX3 *«•**» *»»* •♦*«■•»■*»*»* »•****■*• 
1** i 1.X1M-*,* I'lTIAL STEADY-STATE AT-^SPPEPE ',1«» / TOX 

i tT* ♦ i 'i * * * * ■*'* ♦ ♦ i‘ ♦ * ♦ *r * • » ♦ » > 

Z FCR"AT ( / ^X' VALUEE AT YO* / 5X’P = • ,1 E’ T .3 , * SYnFS/C'Z * 7'^X 
-!*A =*,ni-5,' K-/SEC* 1 SX'TI =*,E1C,2,* •■I***/>«*N =*,tH.I t 
—TSr'X" ="T7 ■■■ 5 X"' ;’;T/Tx *¥Ta .y/S X ' e I-T":?'* ,1 9';'* , T-'€iel ■ 

3-S/ERG’) 

6 FCRYaT t 1 JX’SCLl'IICV OF STEAOY-STATE EDUATIokS’ / 5X 


OTJ- 
7* O 




V ! ' 'a ■ 





tf- JO 


00 

00 

00 


1 •THEK,'fODTKA’«IC VARJAELtS* ! lAK'R’flX 

fR*10X'R-SS*'.CX*P=I)*11X‘P*12X*D* / 5f»(0‘2x'S0L»It MDII-SX 

2»KP’n>.»tCr«<}x'H0N> OIs! 'SX * Off ES/C-«2* SX •|»07C>t7*?X ’6ETA* ) 




FCR^AT t5X*»H£D LARGER DO LOO° LI»=IT FOR » - EXECUTIOA TERrI»;ATER 
1J 


FOB 

?2 format «/ 3X’S3LOTIOV CF STEABY-STATE E-iUATI2Hs’ / 5X 

1 •\OH-0IMEHSrOfAL VAP:ASLES*/'J3X'fi ’nx •VR’YIX’P* riX'DMiX'TMlX'FS 

mnr*' ':a • ?x • k f - • ; ts'^5 — — 

?A FORS^ATU?. =nFE12.A)J 

17 .FORf^AT ( / EX'kACKETIC FIELD* / <3x*R-SO* / tSnP£1i.4}J> 




19 FORHAT <SX*3»<RJ* / ttflPEI *.4j J> 

23 FORBAT f?X*T =*.1PE12.S> 

,'t fCS-AT I / >x“»INinAL E'TERgY tEpsS/RxJ) J — >-'s&Serir"=* .TPEiTrs"; 

1 * , THEF.HAL =* .Eli -S, 

2 *, SPAVITXTIOdAL =*,E1'.5 / SX'IMITIAL PASS 


.;En.i,* &SS/XHZ =*,cn.S,* PARtlCLES/SPm^ 

45 /ORBAT C/3X*REV1SE0 VALUES OF T'"AX A\D R?5*5 */SX 'TPAX =’,1P£9.3 
V CEiiREES* /SX’RPAX =* £9.3,* SOLAR tAMI * /'X * LPX ='.14/ 


■'? FCP’-’AT C/SX'R«DlATI0X*/SX*«3*4X*B0fi D IK * 4» * £96/ CP3-5 * ) 
76 FORMAT C17, 2I1PE12.41) 


RKA 

S3 FORMAT C/* *»* flC COXVERGEXtCE IN STST - ExECOTlOM TcR*«I*i(ATEB ••••) 
95 FORMAT C/5X'ftA0II AT WHICH SOUWO AhD ALFWEK KACH WUP3EKS ARE OWE* 


<24 FCPMAT tl7,3C1PE1 
SET CCHSTAWTS 


OV “ I «3Ct- 

p?i = 1 .67230E-24 
RS = S,957E10 


PlEC2=Pl£C/2. 
RAPCG=. 0174533 


0 Rcvisro A- AX Af.P y: ax, fill X A\0 Y aFP - *S 
XrAX=XfAX+DXT 


X<I«) = CI«-2)»r.XT+TPYH 
_lr {XClK).6E.x7AX-1.E-5 ) 60 TO 1C 

UPITE tIO,13) 

CALL EXIT 


J 

rUI = T«J+(J-l)^DYT 
I F (YC J3 .CP.Y PAX- 1 .E-5 ) 60 TO 12 

"•RITE <10,14) 

CALL EXIT 




iiiiiiilinriirr 


Xf AX = X CLHX-1) 












r«Ax=y (Liy> 

PflMT i5,XF!6X,y«'XX,LPX,Lhy 
L>'>P1 = tJ'X»T 

Aib=PitVrf£7»'ix^t,x}) 

C ylLVES fl yn 

FO = 3.iP'»-'C»Tr 

fto = Sii*i<Pi*<'0/{Fr*D7r7T/i7IT’ 


i. COK>ff = 1, 
GO TO t? 
ts COK«'F=.5 


CO.‘''F=1./.7' 
GO TC fr? 


GO TO c7 
COVt-*F=l./3v'. 

c7 C'd To ri“,"'T5', ■ iJTif"' ~ ' - 

STFA&y-STATc GCLCTICh, ThE":'OlXN*'':C VAHIX’tFS X».C VELOCITIES 
i C"* = SC*FVC’*i 


¥R1TE tlO.iS) 

GO TO 33 
I i = “*i 

CiU. 3CL-? (XC,70,FI.UC,II ,acj 


4-,t 

LL 

5V 1 

‘4 555. 

CS 

1 

556. 

CO 



UTT=US 

S1=OS«»2/t.+1./CPI-T.J-SO 

P^rTE"iT?;3TT 

00 S I=T,1i:'^C 

CC TO C5»,5'3, I5I5I 


EXPlC1*(l./yS-l./YCJ> 


[ 3 01/ » 


faP 10 “I 

**' \ 

! 561, 


60 T5 = YC^{I-l>*ty 

1 

563. 

/>»• 

0Ti=OTl 


,r ^ 

hTT=Os , 


5o4, 

C> 

OS=C . UTT -UTt 


565. 

t ^ 

?■" -'i 


566. 

C6 

Lt^'JS - - - . .. - , 

■' '‘"I- 

56?. 

Co 

FT=f 

• 

563. 

CO 

as='jT-FT7a'-{tj:/y6'»‘3»*tPi-i.>/* T*-oij 




• 


•*>A 







Ln uiui 



j. 

00 

PRINT 

* • 

oc 

im 

CALL ■ 



ps=6s**pr 

TS=P5/0S 

— KT= rrrrrrsT 

P-£S=L 1/AS 

‘1 CO TC a=>,51,i5,52,53,7D,7aJ,I»'*8F 




• 80 TC C50,62),lSTSr 
5l BT<lO=3T»'P$»TrS*»S 

Sf>“Tv> l50,b?i ,I5TST - 

53 eT*0=3T<>»P3*yS«‘1‘' 

CC Tt C50,63J .ISTOT 


n i rt V— -i; 

eo TO C50,62),ISTST 

So C1C.7) l, yS.VS»;,Y58iC,oS.FSa ,DS&,3 T<0 

42 S‘-'ft';=:i*.SS*SOKT<’'TiC»Pl/Z.) 

K^Fi.=Ji/tt£*SCRlC1.-*2./tPI--iTajJJ J 



63 yn,I,1J=TS 

IF ClP?T»<.tr.O) GO TO 73 

“-rF-crsTtTvrc ,t j“t si,=Tc — 

P.ACT\CI> = FADFr CTSC)-uS0'*3 
’3 IF C»S .■E.YF.NO-t.F-SI oO TO 


C ST = A0T-STFTE SCLUTItf., »*A6hETlC F1EL6 
2 X‘'AX''D=T«AX 


•CJ GO TO 7 


KA0TMi=^A5TT<F)*C»ADt. 

PRINT 76,K,RAI»TPN,Fx3TnCIC3 

r’TTrriTrcrF 


DO 22 

rrc~T.TTv:*Tcr-T3 . cxT3>’»i,-T'r;<T 1 72 . " * 
If CXClr.«T.X»AX + 1.F-5> t-0 TO ei 
PTI'jT 2~, XCU^ 

AF.i>l=xcLe»RAPoG 

IF (I»=^X5F,C0.33 i»051=AECT-PIfcC? 

IF (I>'iSF.£r'.53 AFSlsARCT-PlECc 

rt-TI''! Tr;-r .r.Tj jr f- 

j=1,T 

:c TC 

s'‘r?'=irji-. j,T5*»2'" ' ^ " ■ 

bX<1 , J ,-!3--3TMAFGT>/X5C-2*Ce» •** 
sTtn .4 ,T3-£, CC.r(FSC13/Ar.&2*CON'!F 









GO TO 16 

t-> AR62=U ,11**A 

5j:<1 , J = 

TtT'TjTT^^ t t/*CO^ <-^ j 1 > ’* . >/iii=6^ 

cl’ 

aTCI ,J ,i) = CeSCA' G1)*C5.»C<;i<ft‘i6’’**2-3.>/APC^ C N .F 

GC TO 16 _ _ . •■ - 

— ! j ^ 1 ) = ? ^ i * C o A r- i- • 1 1 JTinrWt a R Gi > • i 4757=3>T^ s t ar si ) •^Z*7.SS* 

1COS(AR&1)/t APGtJ 

Id COKTI'^U'^ _ . 

PFI’.T 1 1, t jrci , J,1 J » J=1 fl> 

FILL liKe-PEPcADEHT ARBAtS ' 

24 ms=uc 

CTI^'t "" 

Df. 11 J=1,L'T 

Y r = + CJ-1) *1VT , . ■ ■■■ ^ - i* 

u Cl t J t "t ^ "' 

IF <I3TST.£S.1> £0 TO 35 

«T£=t;Tj_ _ „ .. 

i'Y=us 

if= 2. *yTi-j'’t 

OT=US '’ 

JS=UT-FT/<Ui-<UC/YS**2>«*CFI-1.)/UT’>»FI)^ >,iv- - 

— F='^5ir?r/'Tr*(.i'-/'CT'Y»Tv1i2J 1 /T'-i- ■••'“S' • ^3-.; 

1~ £R -^GS ,S-!jT 5 /t.T .GT,1 »c-i 1 i' "0 

If CAiilF) -LT-1 .£-5) £0 TC °5 

F4. comikuE 
PfI'JT -Z 

CJtt^i^IT 

D{1,JTl> = OOi'<*'5*TU)»*2) 

PCI , jt 1)='' n . j ii > ■*’■^1 _ 

GC' to 15 

j 3 2 Ci,j,i) = cxPCci-a./TCji -i./tc>: 

FCiijfl) “ TCIfJfU _ ^ ■ 

PFDir.C C l , J }=B4r~N<Jl — - 

VfI,J,1> = VC1,J,13 

■jCIfiii-3-*^ 

' i\L’’\Z 1' 3.3tP- 

ztz ,i, 13 = ' ;i.- .1 • 

11 FCifJiii - - 

APG1 =x !■:)■»: 5P;G 







ITtlWACfTfO^I *B61*Ml61-PJECc 
If (IB*6f.eo.7) AlteixAftttt-PtcC2 
00 21 J= 1 ,L**T 


<jC TO Car, 

C OI”CLE FIELD 
57 ASGc=y(J )«»5 


If fl.E«f.1.0R.I.EO.LRJiJ 60 TO 2o 
»23tI-1,JJ=SlKf«?.&1)/Y(J)*«2»CC7.MF 


Ail CI-1 » J } = 

?6 &T(I, J ,1)=2.*C0SCf ss;i)/as62*co.*.»f 
6C TO 21 


APOtt 
1 -25 vjyJ62*TTj) »«* 

“ ' sxritj ,1)=2.«sINtRR61>*C0S(AR£15/A»G2 
IF ti.ES.i. 0 K.l.fci..L' 3 »T~b 0 " TTT 
AidCI-l, J)=SINCARGl)«C0S(»aG1J/Yf J3»«5 
R2Za-1.J) = A22CI-1,JJ*YtJ)*SI\CAKG1J 


■ GO TO 21 

C HEXAPOLE 


EXU, J ,1 J = 3 .*S 1 ' {ARG 13 * ( 5 .-C 0 Si*r 61 )* « 2 - 1 .) /CA.*ARt 2 )*CCN'Xf 

L<Y(I ,J , 1 ) = C 0 SCARG 1 )»( 5 .«COS(X«£l)««' 2 -T. 3 /*BG 2 «CCNPf 


Tf lI.EO.I.GR.I.Ea.LMJ GO TO 21 

AZSlI-1,4> = SIK<*RG1}«t5.«COS<*RGU««2-1 .>/(6.*T CJ) •CONP.f 

A2ZCI-1,j)=AZaCI-1,J)«r(4}*SI«t*a€13 

C OeCAPOLE 

72 ASG 2 =YU )»«7 


BXa,j;i)=25,«cON y F>SIH(AR6 1)«Ct*.«COSIAPG1 J«**/2.-21 .•CCS CARG1 >•* 

12*1.S>/<6.«ARS2) 

EY(I,J ,T3=7.S*CCNPF*CO£tAR61)«C«.»CaS<ARG1 J*»<,/2^5.«COS tARGIJ** 




IF Cl.E0.1.fR.l.E';.LF.X) CO To S’. 

Aifl ( 1 - 1 , J)=CO'<«F» 5 .*{' 53 .»CO':CAi?G 1 )»«A 72 .- 2 l .»C 0 SCARG 13 »« 2 * 1 . 5 >> 


A 2 Za- 1 , 4 >=A 2 Fa- 1 ,J)»Y(J)»SIHtARG 1 J 
21 COUTINUE 

F L :fc ITT ■£ir~Tcr“s ^ 

RADII AT .FlCh .*;5 = 1 A:;f ";A = 1 


00 ?7 I=2,l:*X« 1 
R$ E1=Y0 r- 


ISE 1=1 

IAE 1=1 


•n, 

«». ft * 


IF (R«'2.GE.1.I ISei=2 

ft»-A2=v <1 , 1 , 1 ) • SCFT tC C 1 , 1 , 1 ) / C2X (X ,1 ,1 )««2*ET(1 , 1 ,1 3 **2 ) 1 /AED 
■ If (hrA’;'EtTTr3 rKti='2 
DO Ri j=i,Lvy 
GO TC (•?5.,9G) ,ISE1 

-Tr'RAf5T=T5Tr2 ^ ^ — — — 

aFS2=VCI ,4,1J/SCRTCPCI,J,1)/0CI,J.1) J 
IF (S'.'SS.LT.I,) CO TC S? 






9 • 


/■S.O. 

741 . 

742 . 
— ttjt- 

744 . 

745 . 

747. 
74 3. 
74 5 .-- 

754 

7=1.' 

~? 52 . ■ 

753. 

754 . 


756. 

757. 

~ 75 Tr 


75 y. 

76 n. 

~ 7 FTr 

767 . 
763. 
~7'i4. 
765. 
76 6 .. 
' 7677 - 
'63 . 
76?.’ 

"TTirr- 

771 . 

77 >. 


774 . 

775 . 

ttt: 


77?. 

775. 

~ 77 ^.- 

7 - 


~ 77 T 7 
7c3. 
7=4. 
“Tr?; 
7f 6 . 
737. 


7?i. 

7 = ?. 

"TTr;- 
7? 5. 
756. 


"W 

00 

JC 

dj 

oc 

C3 

-cxr 


Ii£1=Z 

.irt ^ ^ ^.^*1 ^ *077*'C1 . — RJPSl 1 

^1 j sp*'A ^ 

,j, 1 j.s'f, jf. j, 

C-..M A 7 .L 7.1 .■) f.-. 

Ta=i=j 


/cpns 2 -s.‘isi> 


4 -r- ^.lAEl.eS 

4 4 . ' "' ' — 

cl XtI>.«Sp l.SAtl 

'=7 C'if.Ti'U'r t •■• ftt » 


/ C 7 ) 

.2> SC TO 53 


CO 

cc 

TT 

CC 


~w 

cc 


4KSATE Ftti'l FILES 

L-'X--'t^"irCT — 

L-y ^2 = L>'X-£ 

e7y, >-},;rTFF 

75 Cl 2 ,'Zc:,i 7 ;[-- 2 T^rr^ 

no ciP.D 7 cn. 44 i 

’i.3^ > 1 3 Y{ I , j , 1 3 , J=2 ,L 


00 

■ 731 T 

OC 

CC 


~r)Tf f 4 . 4 -:- - -:.t » - » ' i t L 

t i ^ ^ M -. r , L* V *^-U /c V c J >V J - 

. ' XL 

‘C wPI’-c <133 TT^^P. 

^ETnTTTTTl — 


>.t=*-« 2 ,Lvy, 

7 - 4 P- 13 ,J = 1 ,t,. 

L'-rj 


Ptipj, 


13 .Ct* .X.11 

YO.Tva*. 


oc 

?> .-V 


icxcn,i=^,u..x„i 3 ^(yfjj 

?TC , DO , Cl ? ASH 

3 T= Ai 7 ?.'«r 67 TE 


“TF 

00 

00 


STT 

A“: 4 =a:c»* 7 


■ ”X 5 f 3 . .oK" 

l.tMYS.CCAZZ 


. xr' 4 /; 

t J 3 , 1 = 


-3 . 043 . 4.1 
yo'y-sV, — 
1.LrfX»23 ,J 


3. 


>23 ,J=t, 


1 '> 

vO 

~w 


03 

00 


TrAKCY-O.O 

X«AX=XfvAX»ftApD 5 

CXT=CXT*'?ape 3 
T xTF=g 3 Pnx 5 TT^p-r- 

53 4r, X-T,L-X 
XCj 3 =x (X).»! 1 itPDS 
Fo"ii i= 1 ,l‘^^T'~ 


0 3 

'G 


~mr 

03 

S3 


oe IS J=1,L'‘Y 

s ICI,4,13 

V n ., j .TT"= ' 'n'T, j , rj— 

OCX.J,.-) = DCI.J.1) 

_il’'CI,J .2 3 =nXCl.J ,13 



yri, 4 , 2 i= 2 T(:i,T 7 n 

= PCT.J.13 

i.-v5 .-,7 J » f >»»"/cci. j; ir 
• c -J - ■ =3 , J, 1 3 »*i./D Cl , J ,1 3 
"‘=uA = A:-Sl 4 .AfX 2 +A-v- 
:sTrs — 




r r. 


• CPAf DY 3 


~mr 

00 

cc 


2f C 

COMx y., . 

= I'A .1 7 7 =D A I /ChTb^X 

CPAS tT=& XT /I Up ROY 
-ill:.* :^A'‘!CT,ChAPYT 3 r: 




C^A= r; 


G2*.2-4.,ftp5i.*g^2jJ3 J/YtJJ 


T7X ^ im ^ : TT 

10. 

00 

11. 

00 


CALCULATE INITIAL EfcERGIES AMD fASS 
AREAM1=e.f' 

AREAr2=0.C 


aseatt=p.>' 

APEAT2=0.0 


AS£A61=0.0 

Afi£A32=0.Q 

AAEAo— 
AREADI =0.0 
AB£At2=0.0 


DO ?9 J = 2,L^-y 
&T22=3XC2, J ,1)*«2 *£jTC2, 




*REA«1=AaEA"1+.5»<ST21*3T22)*DyT 
AREAT1=AREAT1*.'»fP{2,4-1 ,1) + P(2.J ,1>3«DTT 
AkfcAEI =AREA?rr*.'r» (uRLT+GfiSHJ'Dtr 
AREAa1=AREA0l».5»C0{2.J-1,T)-*0C2,J,lJ3*orir 
ET2l==JT22 


»E1=SR 
DO 2o I=2,I»fXI 

ST21=PX(I,1,13«*2+Eya,1,1>**3 
&REi=o (1 ,1 , i jiTyrrr 
DC 3J J=2,LMy 
eT22=Jxu,j,i)**2+'jy(i, 


AK EAH2 =ASEA?^2+ . 5 • CST2 1 ■»-_<T22 J • C V T 
AREAT2=APEAT2+.^*CPtI, 4-1 .1J + PU , 4 , 1 ) 3*CyT 

AirrToT^T5 s5^r*'T=^iT!PS El -f C A fiXi 3 f “ 

AREAD2=ASEaS 2 + .5* CD C1,4-1 , 13+D 1 1,4,1 ) )*DyT 
3T21=3122 


K-k.S»CA''EAr1 + AAEA**2J »DXT 

7+.5»CAREST1+AeEAT2)*C»T 

feX,^»TAt;£XuTVAA'EA‘125»X.T 
3+.5*(AREADT + A"EA3E)*'»'^ 

«2 
T2 

A!<EAy2=0.0 

AFEAI2=C,2 


AREADI =A RE AD2 
APEAt2=C.O 


AREA?*I 

=A.E 

5»B0»*2»AR£ 


ARE AT I 

=F0» 

f.E**2*A?EAT 

*PI/CPI-1 . 3«1 .P5 

5'=F5'Er 

= «"?■ 

A A 

♦ i> * P \ *?; 5 * f , E '$■/■ V 

APEAC I 

= ARt 

AD»DO*FK«Pi 

**a*1 .55 

AV£AD» 

=AR£ 

AC-I/Pr- 



V.. 











-■ 55^-; CTT“ : 

855. CC 

855. CO 

• FOSTPA»;,C=PLI3,ON=LAXOC1,FL 

SUESOUTIhE LAX<f1 ,!«,DTTK, IL»'X, IFFX, 1 FPYt, 1 FP V2> 

PASAMEIEB t’'.x7 = ;?^,iy7=51 ) 

j-jr 

?5?. :c 

^ ^ 51 . C£ 

STlx’’^ ,yEi<?) , JCIA?, IT?, i) , Ft 1X7, IT?.;. 3 ,5 j IX ■’.IT?, 2),*' ' 

tPClT?, IT7,23 ,L5IS»YCIX73 ,PYfI*r,IT7,2 j,qxtiy7,ir7,C',A2r <1X7,1 Y73 
7 a: 7 aX7,IY73 



IF CKTjiS,M.L 


5*r (#0 

34 G(P, 53=0.0 ^ 


h! 591 . oa 

00 >5 F=1,4 


857. 00 

FN(J',4 3=0.3 


i V 

~35TSr?Tx7=o.-j ■ 


« 654.. C*0 

\= \ 4 


?-9 5. 00 

-Nl = 7 { 

















r 


rOT/r;-bTt>x*T€TTrC>-^ r4 > -o fo t7,U ♦ 

,L53*DTe*tSCS,L)*S<^.,L!l> 

^3t5,L>)/2--3T3X<.«Co<5,L»+rrt6,LJ-G<7,t>-GCt,l>3- 


r 











— 


n = 


.! i 


■U 


t’l- 




TUZ5V 

1026. 

1027. 

-T02T7- 

132<?. 

1030. 

1331. 

1032. 

1033. 

1 C3 ’ S '7~ 
1035. 

1034. 

-Tvrrr 

1033. 

1039. 


“TOCTT' 

1 on . 

1042. 


1047 . 

104r. 


TOroT- 
1050. 
1051. 
~T7JZr 
105 3. 
lost. 

1 055 . 
1053. 
1057. 

Tosrr 

1C5'-. 

1360 . 


TOrSTT 

1062 . 


~nr 

00 

00 


'’0 

CO 

"W 

OQ 

30 


V u 
00 
CO 

"or 

00 

00 


#* 

CO 


CO 

03 


Tf3“ 

CO 

00 


~uU 

oo 

oc 

T3-" 

C3 


To’" 

CC 


(i ,Trr 

0W{3,L) = tst6,LJ*«it 


r 7 j_jV*< 5 ;jyj 

i,LJi72.-BTD**(G<6.l.>-GC5,LJ)-6T0T**tF 

r-in',Tj-i 


— vu VC « y ^ V V ^ c •'*0 i j 9 L J 

If t9,U-Ft2,0-FC3,0)-DT2-CSCt,L>*SC5,L> > 

o.u.L) = cr(47L T7; avLnyT:-sfbMV^re(c~, 

3tT5T-C f(5,C)-FC?,0>+rT2»t'£5,t>»S(2,U ) 


<5.0* 

ci,ijS- 


6 t0(!TI«UE 

7 Ln;L5>F 


"W . _ 

IF (L.E9.1.0R.L.E0.3) T2 =YfJ>»*2 
IF <L.E.S.2> (f2=C<YOF)*V<J 
tf ( L :'F 5'." rin? 2 = f rrrj7*frjKX)vz;3 
YSS=SaRTCT2> 

C».CL.2) = w->(L,?)/4?4(L,1} 


i:* 


6ML.5J » GNtLiiJ/aStL.U 
GH(L.6) * (Cf<(L.6)- 

1 (OfctL.I J« (3K(L,2i»»2*0'«<L,3J*»i2 

F' Cl.,1) = 2J 

£I.<L,1J = Sf^JCfll'C-tC.^I/TSC 


>«. 

-it- 


>*AB02*t 


Tj/rrxrrr 




1043. 

53 


-- 

I 

r 

1044. 

00 

FKCL.Z) * FMIL, 1)*estl.,2>*«»f<L.61/RM02**l&6B2« BTSX2 

' 

t 

1045. 

03 

IF CILOjE.SC.C) CO TO 71 


C0.*.*71*ri,r 

0'«CL,2 5 =(F'' tC.I) ‘4?JCL,5>-S=*C2»C».(L,5 J«0't;<L.t>)/TSi 
f*iCL.31 = G*IIL,2J*VS0 
0«<L,3) = Gf.tL,1J*0«CL,3J* 

1 <Gi:<L,eJ/i»fc02*Jte0a2»t-°TB)t3J>/YS0 


1 " ’^^C3'-’{L,33*QVCL,6>-**i<i_,2?*-N<L,5}5/VSO 

21 = c i *a5{i ' ,*} ' >it.<L ' ;rrF i'eycr.'2)»«g*3'U’L;3 

F»«CL,6> = E1«SS.CL.2J3Ai 302« GX<L,5>*7KtL,5J 
S.VCL.63 = E1«J«<t,33/tSC 

— ♦A-OZ- 'TKXC.E.r'^. ri.', 43 

TF <L.=3.2.tft.L.E7.4> ThZA=XCi) 

IF <Irf£.E9.3.Af.t>.:.E!..L*X' 11 T.-;4 = x<iri 


"TTcT 


IF <lIec. 3) THZ*=t)r<I)*Xtip)>/2 . 


h] 

1163 . 

f 5 

T”ZA = t./TAMIHZiJ 






hi 

roer. 

- 

E-'=c . » jV (L. -';-.,r.rL , : J *T“Z -3 






u 

1 Ct 5 . 

ZL 

SCCL.I J=-:vtL, 1 )'F 2 /Y 5 s. 





t 

6 

1366 . 

-- - 

5 fV(L ,2 J=- 1 '.CL, 13 *<SC./Y 5 .-Z.*. -<C, 2 J** 2 * 1 '«.. 2 >''-Atl-. 3 )»T- 2 F- 




, 

• 


rxrsrz 

—CZ 

ION C L , 3 3 **2 J / YTr-A^Oi* ( isirVrl rj V^-TITsTF 






44 t 

1063 . 

00 

2 THZ*)/YSO 







1067 . 

00 

IF tlLCWa.EI.OJ GD TO 72 




...=« . 


*<■[ * 

T 071 . ■ 

Cl 

ST(L,i; 7 =I\Ct, Tj-ArC?*Cc.»S‘*<L,X>»* 2 * 6 An., 4 ) •JTArC.SJ •iSzF'Cisfc.Al* 







1371 . 


1 CGSC, 5 >-sr.{ 1 . 53 )/EYT 3 /YSC 





» 


1072 . 

r ». . 

72 CUtTiNCE 






M 

•*■ 3 ^ 

T7J7T. 

— 

■■'Si. fc , 21 =-?cTnTT 5 rr 3 T« G\ f ^>Tr*ty CI.T 3 > 7 THi>V** 5 ' 6 ' 2 » CS Ct, 5 > - J 3 .• 






-ri 

1074 . 

CC 

1 Sr.{L, 0 *iC<L. 51 *TJZF 5 /Yti, 





- 


1075 . 

CC 

S‘iJL, 4 >=-S'.<L,i>*TP 2 A 









^-r - — " 






tf- 


1077 . 

IC?"-. 


ToicI 

1051 * 


r 

~T3~ 

*:& 

S3 


If C:F3*%.6=.?) >5 1C 7 
IF fL. •.£.■'> i' 71 It 

46 i" 7 tj 3 - 0?> 3 T *. cTi > : ,"T5 + ysTT^ rr;'j'3'>7£ ■; 


sc TO 7 

:1 If <L.f:E.2> S3 7> 


<w 


77 





















ircvti siGhc=-i. 

CMARTT = SieHC « V { I , Jf !^ l >* CHA«rT 

T -. r^TITt r Mr-£'SS W = tHAlYT ’ 

1 cohHfa *£ 

CC TO 35 


,_J 


NNO £ N=‘j 

39 IF tr.VPR£.i3.0) GO TO GO 

S ' , i ;\* s<s 

J HPRc=C 
C FMSOThl’lG 


If 


HTKSIlT»«T'«SrT^1 
IF C»lfT«)1T.l.T.NTPSK> GO TO 56 

^ ^ ^ ^ ^ ^ 

CALL CAPS r-CItHX) 

56 OTTK = OTT 


CH*PTT=OtT/rHFR&T 

OTT=AJ’i!«1 CC‘-*PXT,CH*.RTT)*ST«»Y 

y’ l ' — — . 

1 A I A V V 

IF <IR3T'<.*.E.2> CO TC 29 


^-IS 


SO It j = i , csa ? 

12 RA0Tl<CCI,4>=5ACF«tPtI.J J,»>1)«COJ**i 


ii'' 15 

15 F.A^T^C (1.JJ=f:ABTfiCC 


CEUai = LMC.'.T *1 
OC 19 J= 1 ,L 5 u =1 

••' F ' ACT ^ir a 3 ' u - • •tj 3 ^ g "& TrtyTOu ? rrx ; a •«- rrrcTi:siiM v 3 %^ 3 T - c £ xcsu " fr , j ,’ s.j 

’ DO ) *•’ * c ^ ■•;• 

* 1 ? ir TO 1 ? 


Kt = 2 

_^Cf» TO 1i_ 

M = 1 
1 i RiTJfe'i 

^ 

5U£!rOUTIV£ l.APS»‘HCIt?*X) 

PARAPETSR CIA7=39,IT7=St) 

— ro ’’ T ' XTT - irriErj '; ym ' 7 Tva ' mnrif 7 ;? 3 ' vvrix 7 ; iTnTr ; TTT *^ 

1P£1*7, I’7,3>,L0IS»nrCIX7),RTa*7,iy7,2)nPXCIX7,IT7,2),«2Ptl)c7,I1f7), 
?SIZC1X7,IV’> 



CCt“*<'CW /CO»t*/ ISMP.Ss'Th.HTPSF.ILOfc'G 
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« = I 1 1 riVI ~ 3 'F t T , 3 , -iT, «■ Cl , ’73 :7 1 " ' ' ' " * 

) = -11,-,) 

^ C'3.F(I,1,\)/MI,3.\3 


0(1 

A=Y2=APC!2'SY(Ir1 ,H>*«2/0tl ,1 ,N3 
FS£2=AKS1*AFX2-<-SlV2 

r- -T? - — - - ■ • - 

£3V'(.5»<ar';+r<;FT(A.-‘S£*F’7»fr-(,>AFG1*'A5Y«:>))3 

£I -*.c = l. 

If 'i'bt'i, 1 ,A’)7rtT’”'5)" iiEAcS^'T. ' 

<;Ka!?AT = 

C3I£\C'>li(I,3,’ >tCfABXT3/Y(1> 

ZJ :TFA7xT7';7T.Cf'F-3A7 fHAf rV^CvAFAT ' ■■' 

Or-.TZM.: 

CFAFfT- S;FT{«5»(AF''2 + StFT(F^f(AFC2*»Z-4.»Ar63"AF'Yi>33> 


IF ('/(I , 1 ,«).L’»T.r ;3 cic?.C=-1. 
CFAF.TT=SIO;,C*V(» ,3,r<3->Cf<FOYT 

I T~ r‘c A A% Y t : * T , c A A t i Y i ^ c ‘t ft f Z r = Cm A Y 3 

ZC:.T3\rc 

CFAF>T=Cf T/ChfFr >. 


T=OYI/CHA?OY 

STTT=AriAl (CPAPXTjCPAAYD^STBY 
IFCOTTT.GT.oTT) go to 3 

r 1 Tj -rn; ^ 3 - : t t , 3 1 - t 7 >« t z * • 

i:r = tif-utt-ittt 


IF (IVLlrt.tC .03 FETL'PJ. 

IF aaEAY.E.-'..F.C?.IB0SY.tQ.73 POIKT 25.IVi.IK 


t<£7>:F\ 

t3C 

-FOOT- Ar ,0=Ft:-i,Z' =0-:hC3.FL 

S'-'eaiJUTiFf 01 0PY( it_«*x,iESo,%,>i1‘5 

: WA'Ch DIKENSIOA Cf ILKY IK cO^fON/COBDOYF 

FASAKETFR (IX7=?9,IY7=51> 

— j ^ f ^ . y 7, i»7,2) VvTI A77lYT,2y7Tn 

'F(I*7. IY7,C},I.D:Ei'Yax?3,er<IX7,lY?,23,EY<IY? ,IV7, 
iACZdA 7.ZY7J 


cci;f;c/ toiSMx.LF'X.L . 

CO/foY/ LKEEP,IL»'Y(IY73.It>(YI 

£': r..r/l3iyxr,irir,XA,-30Y“ 

V -?viT-ri£?t7SAYC?'?S0?ATXj=s sa 
' !»F ZISTUrfPAKCC PRCPAGATFt T« Y^tX 


X - » i Y"^* c 
Z3,SZPtl 


>7* 

XT, IT73 


1CTY?;=0 









VTStsShlT(U(L*I S«l .L.nt )••£»»( L *1 SHI, •?> 
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'fi.Tnr^.L.» 7 i/r(L:::”x....Mi).LT. ffVTjT-sstTr 

rU(L!>ISr*,L,*.>*»J»»fLHS>'l,t.*ll*»2J) >.LT. f**» (V »C J5 

SO T„ 


ueiJHi ■ IDISHI'I 
?6 60 11 L>2.LriEP 
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Lit" » LOIS'tU) 
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1 HlL.LSlS.Hin.LT. » HH . •kfi . « 3 S I ( VTST- So* T ru (L .L SO* . H) tl.LSVi, 

?>!)••;)) ).LT. f "nLfcI$*'TCL»*L 6 IiHr.L >-1 


6 « ?V'L» 5 .L»trP 
LL*LL* 1 

l~-U 

If tU''rlLS>.Ct*. 1 ) fC T 3 ?1 
Lf«.*LDlSH»(LS) 

TT (iFhTTSTrnrrc? ts 2’> 

to 3c I»LIS,Lf»i 

If (*iS( (p(is, I ,M )-p(LS, 1 .Hn/HfLS, I .Hi n .lt. {•p.iHt.issr 

■^iTTnnxT;TV.iJ«*r*» rLsrr.siTi-^ri-sisTruus.i.M-.r.t rts.i .si •• 7 ‘ 

:>).LT. Fo T' ’j 

L5IS *» <LS*n = t 


LIH^Lf H 
2V Ct:«TI7u€ 

If rL 6 i 3 '.*r.LE.i:-'r,tt;:L’».? 4 rii to e 

f.lM 7 

I f ( I ».0) I I\t*2 


ILf « » t 

Lri3*’» * L*'»-1 
t CO Tj (1,J.1<).IL'’»r 

n 7 ;c.rv« “ 

If cLt :» * r ( : >♦; .L t .1 •> uc *o » 

F?:».T I'.i 

XT = =r — ' — — ^ 

ILf»T*< 

CO TC ? 





« CONTltjuE 
GO TG 13 

2 DO ? I=2,LCI5J'J1 

) tc to i 

IF fcO ■^C 3 

XCTY>- = ICTVF-H 


tr-nii 

LBIS«Ya}=LFr«1 
I? {iCTTi*.£e.LK3f-2) 1GF*?T = 3 

rOiNDFIV AT V.'fX 

12 IF ClLFYT.£2*1.PP.I“CrH.cC.»') GO 


.L2IS" 

If ClL»'rCI>.£Q-0) 6C 70 19 

l.S0?1 = tI>^S^yil5-1 

C5-iri=Tni3' Vfiy*1 

l,S<.'=LCII"TC15 

6f T'- C21i2~5f I_IY“ 

■ u < : , L s 't '.T3 = rmf 1 1 , 1 ‘•u o , >1 1 > - 

■ VC I ,LS0fi2,**1I=?.» V<Itt='Ut,K1I- 
0{!,L&'J52,'t1J = ?.* iKI.LSUfc^MJ- 

tiY(l,»£Ui2 ,M) = 2.*3rCl,LiCti,F.1>-i 

F-tr.i-Si;-?,* 1)=?.* p(i,l:.ui:,*.i)- 


ua,LSu*i,«-iy 

VU,t.SU31,«15 

ad.LS’jsij^i 

xiT,tsoSl,M7 

VCI,ISUS1,«.1 ) 
pti.tsur-i.Hi) 


5FfI,LSU 


19 CCHTIKUG 

C %O JNCA»T it XC 

TETTSyu fC-Jf 1 

GO TO (4,',r2flTJ ,I£0X0 
4. DO J=1,LSU3 


V<1»4,N13 =-tfC3.a,NlJ 
Ot1,J,«t) = Of3,J,«13 
/ < 1 , J . •< 1 J =TS T3TJT^1 3 


5Y{1 ,J ,M1 J=-bYC’.J.'11 
9 PCI ,J, M) = PIS, J ,H1 > 


19, 00 

'A - 5 tro'20 J=1,LSUP 

1 0 . 50 

Oa.J.M) =-U{3,J,\1) 

i 1 • 

12. CC 

DC1,J,'D = B{?,J,M) 

i3. ca 

1 bXCI ,J,91)='XC5,J,>«1> 


G = P(3,JrhlJ 

L£2,L?5YtM>=<>.{I 
riT lo IF 
2 GO 23 J=1,LEU'= 

yf1,J,t.13 =-OC3,J,»,13 


BC1, 4,1*1) = P£3,J,K1) 
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'•Oi«.r*»T *T tn» 

17 tt ic it 
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It 

£t;t 
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ii:,i «-»*'•••) ^ 

15 fOSTtT CIS.’C £>J.4)) ■ 
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— »Tr 9 T i-c '•..5 ^ 

<■ f3T'*T **,1Tf10.:,* 
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— 1 ‘Vfr »',rTC.T 7 ," »p/sft’J ■ ~ 

jF»NTxu 

^ 

— T-r x-rrrri'-Tfc ■ ^ 


tit X t.t'i-J 

ITxl 

t'x?./ C-.-I.-SJ 
»*1 . 

TTffiTtC'i ■-■ flTfVT x#*»ic t\Cl SetWTION 

00 t 1 * 1 . ;cc 

6C TC f1t.1t,1«.>,!T 

Tt -.xir-i. 

1« fw.LT.* > CC TC It 
POIM l-.l 


»tTy**r 

1 & S* 4 T** ll‘* ■'/ fV 


l '1 


>1 .tl0**S> 
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APPENDIX B. COHPUTER CODE FOR COHPUTER-OBNERATED PLOTS 

On« outHtandlng advantage of a computing facility aucli cj that 
at tha National Canter for Atmoapheric Research in Boulder » Colorado i 
which is where these sltnulationa were performed, is the capahility 
of quickly producing multi-dimensional plots of the physical variables. 

A listing of the code we used to obtain computer-generated plots 
(some of which are presented in this report) is given here along 
with a description of the variables that must be input to the code. 

The routines in the listing that actually generate the plots are, 

In general, unique to the NCAR Computing Facility; however, the remaining 
logic should be applicable to any facility. 

Since the plots are placed on microfilm, this Introduces the 
poBslbillty of producing computer-generatOC movies. We have produced 
movies of the evolving magnetic field lines for the streamer and 
transient simulations for both 3 * 0.5 and 3 *■ 
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Bescription of Input Variables for Plot Routine PLTFTS for Codes WglBS w4 Mg>Cf (CRAt-l) 


Variable Colums 


Defi nit ion 


Card 1 


Cojawents 


l'nit« 


XDS 


1-A 


Number of data sets to be 
plotted 


NDS should be set equal to the last value of KPTll 
printed out by code MHDHP ' 


I Card 2 | 


ANGH 

1-8 1 

Angular orientation of 
surface plots 

See SCAR 
Vol, 3" 

manual "XCAR Software Support Library, 
for write-up <m Std>iro<i:n; ine SRFACE. 

degrees 

ASGV 

9-161 









Options (Card 2) j 



variable 

Columns 

Definition 

Value 


Conwents 


IPSU 

17-20 

Surface plots of non-di»en- 0 

sional pressure, density, 
velocity magnitude and 
tengjerature 1 

Produce plots - see write-up on SRFACE 
manual. 

No plots 

in above 

IPHT 

21-24 

Half-tone plots of non- 
diraensionsl density, 
teaperatuts; and pressure 

0 

Produce plots - see write-up on HAflOR 
stanuai. 

in above 




1 

No plots 


IPW 

25-28 

Velocity vecl r plot 

0 

Produce plots - see write-up cm VELVEC 

in above 


IPFL 


IPRN 


IPIO 


29-32 


33-36 


37-40 


Contour plots of magnetic 
field lines and the non- 
dimensional density less the 
axabient value 

Half-tone plots of total 
radiation (Cox-Tucker) atni 
the total radiation minus 
the background value 

Plots at t « 0 


1 

0 

1 

0 


0 

1 


manual. 

No plots 

Produce plots - see write-up on OOKSEC In a^bove 
manual. 

No plots 

lr--'iacv plots - see write-up on RA**f>08l in above 
rinual. 

N> p'ots 

Produce i^iots 
No plots 




COHsencs 


ALL contour polar plots. 

Only do calculated field line contour polar plot. 

Include on plots. 

Not on plots. 

No dock. 

Indude dock. 

No flash buffers tised. 

Use flash buffers. Set all variables to 1 except 
lPFIXKIPTO-IPRSA-0. 

less than 1000 franes so itse usual DISPOSE card 
procedure. 

More than 1000 fraaes so reaove DISPOSE cards and 
use Fortran call to dispose. 
















Listing of Computer Code PLTPXS 


SELT.IOL .S*M 
ELT 8B1 S7S0T 

1 . 

2. 

J. 

Am 

5. 

6 . 

i.. „ , 7 . 

8 . 

9. 



14. 

00 

V ■ 

15. 

00 

00 


17. 

00 

w 

18. 

00 


19.,,., „„ .... 

.„.00 


20. 

00 

V - . 

21. 

00 


77 . 

0.0. 


23. 

GO 

V 

24. 

00 


. .. 25. 

00 


26. 

00 

^ • 

27. 

00 

T." 

. 28. 

^00. 


29. 

00 

W -- 

30. 

00 


PLE, .STOP 

C 04/25/80 11:18i4A C->0> 

00 •JOB. 5346, 35881000, STEINOLFSON 

00 .FORTR»N,D=PLIB,Oh=PLTCT2,AF,FL 

00 P«0G«*« PLTPTS 

00 C 1. I t 4M D JJ BUS T E E AT LEAST .*?_ L»«6E AS J»*f*FS.aEA6 , 0W_ TAP|j_ 

00 C II.SE.LP-X-2, JJ.6E.LHT, 1I«IX7. JJ*IT7 

00 C 2. DIMENSIONS OF 2 BUST BE EOUAl. TO ACTUAL HtllHBEA OF POINTS PLOTTED, 

00 C -B7 AND N7, R7=IRXP, N7*JBXP 

00 C 3. DIBENSlONS OF NORK BUST BE AT LEAST AS LARGE AS 2 «b7.N7*I|7*N7 

00 C A. lHr*2*IBXP»JB*P*IBKP^JB*P, 1IV*III«P/2*1 , Jj V*JBXP/2*1 

00 C 5. N »1FL«3«JHXP 

00 PARABETER {Ix7=37,TY7*41 ,IBKP*3?, JflXP*At , t«*2ATA,lIV»t9, JJV»21 , 

00 mBFL*123) 

00 „ _ _ real LABI, lab . _ . 

00 DIRENSION P (1X7, ir7J,DCIX7, 1771, V<IX7,iy7J,UCIX7,iy7), axt 1X7, IT7), 

00 IBY tIX7,lY7) ,X(1X7) ,LTtIX7), 7(1771, FLtIX7,IY7J, 

ao 21AC.U2J..DAClT7.1..aADJNClT7> 

00 DIHENSION Z(IBxP,JKXP),M0RX(IHK),VU(lIV,Jjyi,VV(IIV,JJV),SPV(2), 

00 lLA91(31,LAB(6J,2CU«XP,IBXP»,2Z(JPXP,IBXli'3,WUC(JJ»,llV), 

,„.00 2WVC<JJV,lIV),XXLaX7),YTL(lX71,XXOax7),Tm'Uax7),VUD(JJV,llV), 

00 3VV0rjJV,IIVa,XXB(NBFL),7YB(NBFL) 

00 COBNON /TRANS/T1,Y2,X1,X2 

0 0 56 FORBA T (■ . *• INCREASE NBF L ••• ’1 

GO A FORBAT (1-(/2F3.0,6I 41 

00 3 FORBAT (’IINHUT -*//3X*N0S **IA / 

00 _ _ , 23X*ANSH .= ’E11.A / SX'ANGV s'EH.A / 3X*1PSU **,I3 /.SX’IPHT «*, 

00 3I3/3X'IPVV **,I3/3X'1PFL ** ,I3/3X* IPRN *• ,I3/3X* IPTO *',I3J 

00 31 FORBAT ('TIME **,1PE10.3,* BIN., *1 

^00 SPjfi 1 13 0 .i)_. ^ ^ 

00 SPV(21=0.0 

00 IPXP**IBXP/2 

00 IF (2«IRXPV.LT.IBXPl IBXPV»IBXPV*1 

00 JBXPW=JBXP/2 

00 IF (2»J«XPV.LT.JRXP) JBXPV*JBXPV*1 

0 C=a-. _C PLOT OPTIO N - W ATCH Z DIBENSION 

00 1PL0T=1 



. . N0ST=1 . . . . 

P1EC=3. 1415926 
RLTIC=.05 

RA.PDS=;.,01Z4 533 - 

READ (101 71BE,B,N, C (P (I , J 1 , D (1 , Jl ,U(I , J 1 , V (I, J> ,3X( I, J 1 ,BY(I , J 1 , 
1I=1,S1,J*1,N1,XRAX, YO,YBAX, (X(I1 ,1*1 ,Bl , (7 tjl , J=1 ,N1 , {(FLdmJl » 

2I=1,Bl,Ji1,N},TO,DO,CRAON , . ^ 

Tl=Y(11 

T2*T(JBXP}-T1 

AB»rltT2 

DO 49 I=1,IBXP 
X(I1=XCI1«RAP0G 

.XXL(I1=Y(1)«C0S{PI£C-X(I11 ^ 

YYL(11=YC11*SIN(PIEC-X(I)1 

XXUtIl=YCjBXPl*COSCPIEC-X(Ill 

_4AjTUdl=^^t(JJEXEi?SINC£lECj=XOJll__ 

y1T2*Y1*Y2 
X2=X(IRXP1-X(11 
_ X1*0.0 _ 


V 


J JJ J * 





56. 

57. 

. . 58. 

59. 

60. 

61 . 

62. 

63. 

..64 , . _ 

65. 

66. 

67. 

OS 

00 

OO . 

CO 

00 

00 

00 

00 

00 .. 

DO 

00 

00 


6B » 

00 


69. 

00 


70. 

OO 


71, 

30 


72. 

00 


73 

00 


74. 

00 

■m 

75. 

00 


76. . 

00 


77. 

00 


7f . 

00 



79* 

.00 


80. 

00 

m- 

81. 

00 


82* . ■ 

. . 00 


63. 

00 

9 

84. 

00 

.. . 

35. 

„oo 


86. 

00 

• 

87. 

00 


88. 

. . 00 ...... 

* 

89. 

00 

• 

90. 

00 


91 • . 

00 


92. 

00 


93. 

00 


94. . . 

... .00 . . 

, r 

95. 
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00 ... 


110. 
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112. 

33 


read C5f43 «OS,Af<GH,ANGV,IPSU,IPHT,lPVV,IF'FL,IPR«»IPT» 

PRINT 3, NDS,AN6H,ANGV,lPSU,lPHT,IPVV,IPfL,IPRN,lPT0 
- -DO 22 J=1 ,K . 

TACJ)=P(1,J)/0a,J3 
22 PACJ)=0C1,J) 

IE. aPRN.E<i.U GO T0„A4 _ ... 

l>0 25 J=1 ,t« 

25 RADTNt J)=RA0fNCTA{Jj*T05*<0AO)*00J«*2«CRA0N 
. 44 IF CIPTO.EC.I) GO TO 16 „ ... 

ENCODE 124, 31, LABI) TIME 
LAB(1)=LAB1 CD 

XAaI2X=LAfiJ.C25 ; 

LAaC3)=LAB1 (3) 

,G0 TO 26 
46 iPTO=1 

16 read CIO) Time, lx,cltcj),j=t,lx> 

ENCODE C24,31,LAS1) TIME 

...™.:._LABa).=LAal CD .. _ 

LAEC2)=LAB1 C2) 
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__.DC 14. 1=1. LX . 

L =LTCI) 

READ CIO) CPa,J),DCl,J),HXCl,J),EVCI,J),aCI,J),VCI,J),FLCI, 

~~14 CONTINUE " 
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.IF CIPSU.EC.1) 60 TO 13 _ 
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__.17_IF_CfiDST^£C.1) 60 -TO 7 . . . .. . . 
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GO TO 6 
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6 CONTINUE 
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LABC5)=1H 
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60 TO 36 
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GO TO 36 
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LA5C5)=3HURE 

LABC6) = 1H . ... . . ... 
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CALL EZSRFCCZ,I«4XP,IKXP,ANGH,ANGV,kORK) 

7 continue 
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: 


1*3. 

00 


60 TO 28 



1**. 

00 


30 Z(I,J)»ZCI, J)-RAOTN(J) 



_ „ 1*5.- 

00_ 



.28 CONTINUE, . . 
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'j 
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256. 

00 

74 

continue 


257. 

00 

88 

IPFL1*1 

# ^ 

258. 

00 

73 
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CALL LINE (0.0,0.0,AB,0.03 
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90 

CONTINUE 

# 
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266. 

00 
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,,4'; 
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00 

61 

CONTINUE 


270. 

00 

62 

CALL CURVECXXL.TTL.inxP) 

... ., .^ 

271. 
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CALL CURVE (XXU, YTU.IPXP) , . 
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272. 

00 


GO TO (71 ,72) ,IPFL1 
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00 

71 
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332. 

00 

• ^ 

333. 

00 
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335. 

GO 

336. 
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338. 
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• 

339. 
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340, 

00 


T1=XCL) 

XX«C1 J=ai*C0StPIEC-T1 3 
irY»iCl J=R1»SINCPIEC-T1 > 

1CO!«P=0 

BR1-3y{t,T) 

- ~.SXT=BXC1.,1J 

CO 51 K=E,NKfL ■ 

KK=K 

- br=dl 

IF lABStgTI ),6T.1 .E-5’AfiS(8Rl> J &R=31.»AES<;BR1 /BT1) /SiRKI .-FCBRl / 
1BT15**2} 

XF JCBRa.tX.JB.03..!>Sr-CR 

01=CL/R1 

IF CABS tBRI J .GT.1 . E-5*ABS CBT1 >3 DT= ABS CDR-« BT1 /BR1 3 /R1 
- IF CBT1,LT.0.0> BT=-CT 
S1=R1+D8 
T1=T1+CT 

IF_CB1.ST_.ra3>, 60 ..T 0 63 

R1=YC1> • 

ICOHP=1 

63 IF CR1 .LT.YCJMXP33 GO TO 6tf 
R1=YCJMXP3 
ICOHP=1 

&4_Xf_CTl.ET..JC.a>i._6a TC3 65 

T1 = XC13 ~ ' ■ 

1C0PP=1 

65 IF CTI.LT.Xtl«XP>>GO TO 66 
T1=XCIMXP3 
1C0MP=1 

— 66 JCXRCK3=jat.c.OSCei£C-T13 _ 

yTHCK3=R1*SIN(PIEC-T13 ” ^ ' 

IF CIC01P.t0.13 60 TO 57 

DO 52 I=2,HXP 

II=I 

IF CXCI3.GT.T13 GO TO 53 
_52_rONTIKUE„„ 

53 1=11-1 ~ *“ - — ' ' 

IP1=I-t-1 

00 54 J=2,JKXP 

>J = J 

IF CYCJ3.GT.H13 GO TO 55 
— 54_ COMlfJUE _ „ 

55 4=J4-1 “ ' ' ' ' 

JP1=J+1 

. , 8P1=3Ya,J3FCP1-YCj3)*<BYCI,JP13-BYCl,JJJ/{YCJP1J-TC4)J + CTl-xa>3* 

1CBYCIP1,J3-8YCI,4>+Cg1-YCJJ)*CBYCIP1,JP1 J~BYaPt,43-BYCl,jp1)+ 
23YC1,433/CyCjP13-YCJ333/CxC1P1J-XCI3J 

BT1=BX(X.4> + {H1_-YC4>3»C9X.CI ,4P1 J-BX<I , J> 3/ C YC4P13-YC J J l-^CTI-Xf I) 3 * 

1 CEXCIP1 ,JJ-BX<I,J)+CR1-YCJ33*CBXCIP1 , JP13-BXCIP1 ,43-BX(I,JP13V 
2BXCI,JJ)/ttCJP13-y(j)53/CXClPl)-xU33 
-_51 COATINUE 
PPINT 56 

57 call CL'PVE CXXM,tVK,KK3 
_ 70_ COMTinuE, 

CALL FRAKE ~ ' " 

IF ClPTO.cC.03 GO TO 46 
21 liCST=NDST*1 
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341. 

342. 
^ 343 .. 

344. 

345. 
344.. 

347. 
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349.. 
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351. 

352.- 

353. 

354. 
„ 355.. 

356. 

357. 
358... 

359. 

360. 

-_361.. 

362. 

363. 

364._ 

365. 

366. 
._367._ 

368. 

369. 
370.. 

371. 

372. 
-_373.. 


00 IF <N0ST.LE.N0S) 60 TO 16 

00 END 

.00 fUNCTlOM »«»FNCTS) ... „ . 

00 OIMENSION TC4«),a<49) 

00 04T4 T/9.2e3.1.Eo.1.0SE4,1.1Se4,1.26c4,1.«SE4.1.BS£4,2.E4.2.2E4, 

JB.0 12.3E4..Z.65 EA.l, £.4 ,i. 5.E4 ,.4_. J.E4.4..9.Eit,5»«4, 6.. J 7f.4 ,.7..0;E4 .0.5E4. 

00 21.1E5.1.3E5«1.65ES,1.9E5,2.17CS,2.63CS,2.8SeS,2.98E5.3>3ES* 

00 33.9E5,4.7E5,5.63eS,6.2E5.6.9ES*7.4E5,«.t>7E5.9.08ES,1.C6E6,1.3E6, 

JO 41.75E6.2..2BE6,3..ie6,4..3E6,7.03E6,9.85E6,1.25E7j2.23E7,3.7E7, 

00 56.02E7,1.E8/ 

00 DATA Q/1.1E-24,3.6E-24,8.0SE-24,1.65E-23,3.5E-23,5.85e-23. 

JO U.35tr2ia.JE=22J.07E-;22,9.5E.-:23,8.«e-23,1.E-22,t.2E':2.2t1j.0]lC-7.?. 

00 22.23E-22,3.1SE>22.4.3E-22,6.05E>22.7.25E-22,8.04e-22.8.0SE-22. 

00 38.3E>22,9.3E>22.1.02E-21>1.E-21.8.eE-22,7.6E>22.6.02E>22,S.1E-22, 

00 ^il4.,7fr22i4.2E-22,3.8E-22,2.84E-22,1.9E-22,1-24E-22,9.3E-23» . . 

00 57.1E-23,6.E-23,5.7E-23.5.3E-23,4.9E-25,3.25E-25,2.4E-23,1.9E-23, 

00 61.7E-23,1.9E-23,2.25E-23,2.9e-23,3.47e-23/ 

JOl 3 FOBWAT (2X'«»« T.GT.1.E 8. AA D I ATI 0|l_ EATAAFIQ LaTEb «»•■) 

00 IF CTS.6E.Tf1)> 60 TO 1 

00 RADFN«0.0 

JO AETUIlll _ _ ... ... . 

00 1 IF (TS.LT.TC49)) 50 TO 2 

00 PRINT 3 

Jfl RADE N» QC47>*C 8C49)-a(47»» CTS-Tt47 >>/ tTC49)-T{ 47) > 

00 RETURN ■ 

00 2 DO 4 I«1,43 

.00 IP=I*1. _ . 

00 IF CTS.SE.TCD.AND.TS.UT.TCIP)) 60 TO 5 

00 4 CONTINUE 

JO 5 8ABFN=aCI)*(aClP)-a (l))«CTS-TCI>>/CT< l P>»T U>> 

00 RETURN 

00 END 



END ECT. ERRORS: NONE. TINE: 3.532 SEC. INa6e COUNT: 373 
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APPENDIX C. SAMPLE RUN 

A llaclng of a typical eonputec run for a coronal atraamer 
aimulatlon la praaented. lot Input varlablea are Hated first 
followed by the initial ambient atmoapherc and the time-dependent 
solution. The time-dependent solution la given at selected time 
Intervals along radial lines at selected meridional Increments. 
For this example the solution Id only shown for the first time 
increment. 
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